General Solution to a Simultaneous Linear Equations

Description

General Solution to a Simultaneous Linear Equations

Usage

LinEq(A, c, Ginv = 0)

Arguments

Details

The general solution, x0, to the equation

    A %**% x = c
 

is given as

    x0 = ginv(A) %*% c  +  nullA %*% y
 where
    nullA = diag(ncol(A))-ginv(A) %*% A
    and
    y is any conformable vector.
 Note that, since A %*% nullA = 0 and  c = A %*% z, z is any,
    A %*% x0
      = A %*% ( ginv(A) %*% c + nullA %*% y ) = A %*% ginv(A) %*% c
      = A %*% ginv(A) %*%   A %*% z = A %*% z = c .
 

When A is square, ginv(A) will be calculated by matSwp function by default.

Check mad and max_nullA:
mad >> 0 if the system is inconsistent.
max_nullA >> 0 if the solution is not unique.

Value

A list of x0, nullA, mad, max_nullA
where
mad=max(abs(A%*%x0-c)) for checking inconsistency, and
max_nullA=max(abs(nullA)) for checking uniqueness.

Examples

# A %*% x = c with full rank A: x = inv(A) %*%c is the unique solution
A=matrix(c(
 1,1,-1,
 4,3,-1,
 6,5,-4
), 3,3, byrow=1)
c=A%*%c(1, 3, 0)

# unique solution
LinEq( A, c )


# rectangular A
# drop the last column: too many equations.
A1=A[,-3]

# general solution:  nullA is Zero (unique solution).
sol=LinEq( A1, c )

# solution 1
set.seed(1703)
y=sample(0:9,ncol(A1),replace=1)
Print(sol$x0, sol$nullA, sol$max_nullA)
x=sol$x0 + sol$nullA%*%y
Print( A1, x, A1%*%x, c, y, sol$mad )

# solution 2 == solution 1
y=sample(0:9,ncol(A1),replace=1)
x=sol$x0 + sol$nullA%*%y
Print( A1, x, A1%*%x, c, y, sol$mad )


# drop the last row: too many unknowns (no unique solutions).
A2=A[-3,]
c=c[-3]

# general solution
sol=LinEq( A2, c )

# solution 1
set.seed(1703)
y=sample(0:9,ncol(A2),replace=1)
Print(sol$x0, sol$nullA, sol$max_nullA, sol$mad)
x=sol$x0 + sol$nullA%*%y
Print( A2, x, A2%*%x, c, y )

# solution 2
y=sample(0:9,ncol(A2),replace=1)
x=sol$x0 + sol$nullA%*%y
Print( A2, x, A2%*%x, c, y )



# non full-rank A with c=0: Homogeneous System of Linear Equations
# (no unique solutions)
A=matrix(c(
 1, 1,-1,
 4, 3,-1,
 1, 1,-1
), 3,3, byrow=1)
c=rep(0,3)

# general solution
sol=LinEq( A, c )

# solution 1
set.seed(1704)
y=sample(0:9,ncol(A),replace=1)
Print(sol$x0, sol$nullA, sol$max_nullA)
x=sol$x0 + sol$nullA%*%y
Print( A, x, A%*%x, c, y, sol$mad )

# solution 2
y=sample(0:9,ncol(A),replace=1)
x=sol$x0 + sol$nullA%*%y
Print( A, x, A%*%x, c, y, sol$mad )


# Warning: x0 below is NOT a solution because the equation is inconsistent.
A=matrix(c(
 1,1,-1,
 4,3,-1,
 0,0, 0
), 3,3, byrow=1)
c=c(1, 2, 3)
sol=LinEq( A, c )
Print( A, sol$x0, A%*%sol$x0, c, sol$max_nullA, sol$mad )