How to solve linear system $A x= b$ when matrix $A$ is non-invertible?

1.2k Views Asked by At

I need to solve an equation of the form $Ax=b$. Since $A$ is a singular matrix so it doesn't have an inverse, is there a way I can calculate $x$ by solving this equation?

2

There are 2 best solutions below

0
On

if $A$ is singular, you either want to solve one of these two problems: \begin{align*} & \label{pb1}\tag{1} \min_x\|Ax-b\|^2 \\ & \label{pb2}\tag{2} \min_x\|x\|^2:Ax-b=0 \end{align*} For the problem \eqref{pb1}, you want to use Newtons' method for optimization. For the problem \eqref{pb2}, you want to use Lagrange multipliers.

2
On

You need to reduce the augmented matrix $[A, b]$ into RREF (reduced row echelon form). This will simplify finding the general solution. Below are procedures written in Visual Basic to do just that.

Public Sub solve_rd_system(q, b, m, n, mt, x0, xa, ierr)
Dim lones(200) As Integer
Dim jt(200) As Integer
Dim i, j As Integer


' Input matrices q and b of the system q x = b where q is m x n
' Solution vector is of the form:    y = x0 + xa t
' where x0 is m x 1 , xa is (m x (mt) ), and t is arbitrary vector of 
' dimension ((mt) x 1), the first (mt) columns of matrix (xa) span the null
' space of matrix q, i.e. they form a basis of the solution of q x = 0  

eps = 1.E-8
ierr = 0

For i = 1 To m
q(i, n + 1) = b(i)
Next i


Call reduce_matrix(q, m, n + 1, det)

   For i = 1 To m

   lones(i) = -1
   For j = 1 To n
      If Abs(q(i, j) - 1) < eps Then
         lones(i) = j
         Exit For
      End If
   Next j

   Next i

' solve for x0vec

For i = 1 To n
x0(i) = 0
Next i

For i = 1 To m

j1 = lones(i)
If j1 <> -1 Then
   x0(j1) = q(i, n + 1)
Else
   If Abs(q(i, n + 1)) > eps Then
      ierr = 1
      Exit Sub
   End If
End If

Next i
tcount = 0
For i = 1 To n

   For j = 1 To m
       If i = lones(j) Then GoTo lnext_i1
   Next j
   
   tcount = tcount + 1
   jt(tcount) = i

lnext_i1:
Next i

For i = 1 To n
For j = 1 To tcount
xa(i, j) = 0
Next j
Next i


For jtndx = 1 To tcount
xa(jt(jtndx), jtndx) = 1

For i = 1 To m
j1 = lones(i)

If j1 <> -1 Then xa(j1, jtndx) = -q(i, jt(jtndx))

Next i

Next jtndx

mt = tcount

End Sub

Public Sub reduce_matrix(a, ByVal m As Integer, ByVal n As Integer,
  ByRef det As Double)

Dim i, j, k, k1, krow As Integer
Dim amax As Double, kmax As Integer

eps = 1.E-8
krow = 1
det = 1
For i = 1 To n
    
k1 = 0

amax = Abs(a(krow, i))
kmax = 0
For k = krow To m

   If Abs(a(k, i)) > eps Then
      If Abs(a(k, i)) >= amax Then
         kmax = k
         amax = Abs(a(k, i))
      End If
   End If

Next k

If kmax <> 0 Then
   k1 = kmax
End If


If k1 = 0 Then
   det = 0
   GoTo lnext_i
End If

If k1 <> krow Then

' interchange rows k1 and krow

   For j = 1 To n

      t = a(k1, j)
      a(k1, j) = a(krow, j)
      a(krow, j) = t

   Next j
   det = -det
   
End If


'  normalize krow-th row

   t = a(krow, i)

   If Abs(t - 1) > eps Then
   
   For j = i To n   
   a(krow, j) = a(krow, j) / t   
   Next j
   
   det = det * t

   End If
   
' eleminate elements in column i

For i1 = 1 To m

If i1 <> krow Then

   factor = -a(i1, i)
   If factor = 0 Then GoTo lnext_i1
      
   For j = i To n   
   a(i1, j) = a(i1, j) + factor * a(krow, j)   
   Next j
   
End If

lnext_i1:
Next i1


krow = krow + 1

If krow > m Then Exit For

lnext_i:

Next i

End Sub