System of Linear Equations

System of Linear Equations

\[ {\bf A x= b} \]

Gaussian elimination

### Setting up library
library(matlib)

## Example 41
A <- matrix(c(1,2,1,1), nrow = 2, ncol = 2)
A

     [,1] [,2]
[1,]    1    1
[2,]    2    1

b <- c(4, 5)
b

[1] 4 5

Solve(A, b)

x1    =  1 
  x2  =  3 

solve(A, b)

[1] 1 3

plotEqn(A,b)

  x[1] + x[2]  =  4 
2*x[1] + x[2]  =  5 

A <- matrix(c(1,-2,1,-2), nrow = 2, ncol = 2)
A

     [,1] [,2]
[1,]    1    1
[2,]   -2   -2

b <- c(4, 5)
b

[1] 4 5

Solve(A, b)

x1 + x2  =  -2.5 
      0  =   6.5 

solve(A, b)

Error in solve.default(A, b): Lapack routine dgesv: system is exactly singular: U[2,2] = 0

plotEqn(A,b)

   x[1]   + x[2]  =  4 
-2*x[1] - 2*x[2]  =  5 

A <- matrix(c(1,-1,1,-1), nrow = 2, ncol = 2)
A

     [,1] [,2]
[1,]    1    1
[2,]   -1   -1

b <- c(4, -4)
b

[1]  4 -4

Solve(A, b)

x1 + x2  =  4 
      0  =  0 

solve(A, b)

Error in solve.default(A, b): Lapack routine dgesv: system is exactly singular: U[2,2] = 0

plotEqn(A,b)

 x[1]   + x[2]  =   4 
-x[1] - 1*x[2]  =  -4 

## Example 42
A <- matrix(c(1,-2,-1,2,3,2,3,-2,1), nrow = 3, ncol = 3)
A

     [,1] [,2] [,3]
[1,]    1    2    3
[2,]   -2    3   -2
[3,]   -1    2    1

b <- c(6, -1, 2)
b

[1]  6 -1  2

Solve(A, b)

x1      =  1 
  x2    =  1 
    x3  =  1 

solve(A, b)

[1] 1 1 1

plotEqn3d(A,b, xlim=c(0,4), ylim=c(0,4))

Data example

One of the most classical problems in optimization is the transportation problem.

Let \(x_{ij}\) denote the amount shipped from a production site \(i\) to the distributor \(j\).

## Example 38
M <- matrix(c(1,0,1,0,1,0,0,1,0,1,1,0,0,1,0,1,475,489,542,422), nrow=4,ncol=5)
M

     [,1] [,2] [,3] [,4] [,5]
[1,]    1    1    0    0  475
[2,]    0    0    1    1  489
[3,]    1    0    1    0  542
[4,]    0    1    0    1  422

M2 <- M[,-5] 
M2

     [,1] [,2] [,3] [,4]
[1,]    1    1    0    0
[2,]    0    0    1    1
[3,]    1    0    1    0
[4,]    0    1    0    1

## Example 43
A <- matrix(c(1,0,1,0,1,0,0,1,0,1,1,0,0,1,0,1),nrow=4,ncol=4)
b <- c(475, 489, 542,422)
Solve(A, b)

x1     - 1*x4  =   53 
  x2     + x4  =  422 
    x3   + x4  =  489 
            0  =    0 

Data on supply of Malta Guinness

We are trying to allocate transportation of goods (supplies) to other locations according to their demand.

### Practical Applications
A <- matrix(c(1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 0 ,
0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 ,
0 , 0 , 0 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 1 , 0 ,
0 , 0 , 0 , 0 , 0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 , 0 , 0 ,
0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 1 , 0 ,
0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 , 0 ,
0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 1 ,
0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 , 0 ,
0 , 0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 1 ,
0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 ,
0 , 0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 ,
1 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 ,
0 , 0 , 0 , 0 , 0 , 1 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 ,
1 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 0 , 1 , 0 , 0 , 0 ,
0 , 0 , 0 , 0 , 0 , 1), nrow = 11, ncol = 18, byrow = TRUE)
 b <- c(1298, 1948, 465, 605, 451,  338,  260,  183,  282, 127, 535)
Solve(A, b)

x1                    - 1*x11 - 1*x12 - 1*x13 - 1*x14 - 1*x15 - 1*x16 - 1*x17 - 1*x18  =  -1483 
  x2                    + x11                                                          =    605 
    x3                          + x12                                                  =    451 
      x4                                + x13                                          =    338 
        x5                                      + x14                                  =    260 
          x6                                            + x15                          =    183 
            x7                                                  + x16                  =    282 
              x8                                                        + x17          =    127 
                x9                                                              + x18  =    535 
                  x10   + x11   + x12   + x13   + x14   + x15   + x16   + x17   + x18  =   1948 
                                                                                    0  =      0 

Infinitely many solutions!

Echelon Form

### Setting up library
library(matlib)

## Example 44
A <- matrix(c(0, -1, 1, 0, 1, 1, 0, 1, 3, -4, 2, 0, -1, 0, 4, -4), 4, 4)
b <- c(1, 1, 5, -2)
showEqn(A, b)

 0*x1 + 1*x2 + 3*x3 - 1*x4  =   1 
-1*x1 + 1*x2 - 4*x3 + 0*x4  =   1 
 1*x1 + 0*x2 + 2*x3 + 4*x4  =   5 
 0*x1 + 1*x2 + 0*x3 - 4*x4  =  -2 

echelon(A, b, verbose=TRUE, fractions=TRUE)

Initial matrix:
     [,1] [,2] [,3] [,4] [,5]
[1,]  0    1    3   -1    1  
[2,] -1    1   -4    0    1  
[3,]  1    0    2    4    5  
[4,]  0    1    0   -4   -2  

row: 1 

 exchange rows 1 and 2 
     [,1] [,2] [,3] [,4] [,5]
[1,] -1    1   -4    0    1  
[2,]  0    1    3   -1    1  
[3,]  1    0    2    4    5  
[4,]  0    1    0   -4   -2  

 multiply row 1 by -1 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1   -1    4    0   -1  
[2,]  0    1    3   -1    1  
[3,]  1    0    2    4    5  
[4,]  0    1    0   -4   -2  

 subtract row 1 from row 3 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1   -1    4    0   -1  
[2,]  0    1    3   -1    1  
[3,]  0    1   -2    4    6  
[4,]  0    1    0   -4   -2  

row: 2 

 multiply row 2 by 1 and add to row 1 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1    0    7   -1    0  
[2,]  0    1    3   -1    1  
[3,]  0    1   -2    4    6  
[4,]  0    1    0   -4   -2  

 subtract row 2 from row 3 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1    0    7   -1    0  
[2,]  0    1    3   -1    1  
[3,]  0    0   -5    5    5  
[4,]  0    1    0   -4   -2  

 subtract row 2 from row 4 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1    0    7   -1    0  
[2,]  0    1    3   -1    1  
[3,]  0    0   -5    5    5  
[4,]  0    0   -3   -3   -3  

row: 3 

 multiply row 3 by -1/5 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1    0    7   -1    0  
[2,]  0    1    3   -1    1  
[3,]  0    0    1   -1   -1  
[4,]  0    0   -3   -3   -3  

 multiply row 3 by 7 and subtract from row 1 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1    0    0    6    7  
[2,]  0    1    3   -1    1  
[3,]  0    0    1   -1   -1  
[4,]  0    0   -3   -3   -3  

 multiply row 3 by 3 and subtract from row 2 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1    0    0    6    7  
[2,]  0    1    0    2    4  
[3,]  0    0    1   -1   -1  
[4,]  0    0   -3   -3   -3  

 multiply row 3 by 3 and add to row 4 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1    0    0    6    7  
[2,]  0    1    0    2    4  
[3,]  0    0    1   -1   -1  
[4,]  0    0    0   -6   -6  

row: 4 

 multiply row 4 by -1/6 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1    0    0    6    7  
[2,]  0    1    0    2    4  
[3,]  0    0    1   -1   -1  
[4,]  0    0    0    1    1  

 multiply row 4 by 6 and subtract from row 1 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1    0    0    0    1  
[2,]  0    1    0    2    4  
[3,]  0    0    1   -1   -1  
[4,]  0    0    0    1    1  

 multiply row 4 by 2 and subtract from row 2 
     [,1] [,2] [,3] [,4] [,5]
[1,]  1    0    0    0    1  
[2,]  0    1    0    0    2  
[3,]  0    0    1   -1   -1  
[4,]  0    0    0    1    1  

 multiply row 4 by 1 and add to row 3 
     [,1] [,2] [,3] [,4] [,5]
[1,] 1    0    0    0    1   
[2,] 0    1    0    0    2   
[3,] 0    0    1    0    0   
[4,] 0    0    0    1    1   

## Example 45
A <- matrix(c(0, -9, 4, 4, 5, 6, 3, 8, 4, 6, 0, 4, 7, 7, 5, 12), 4, 4)
b <- c(-9, -17, -3, 0)
showEqn(A, b)

 0*x1 + 5*x2 + 4*x3  + 7*x4  =   -9 
-9*x1 + 6*x2 + 6*x3  + 7*x4  =  -17 
 4*x1 + 3*x2 + 0*x3  + 5*x4  =   -3 
 4*x1 + 8*x2 + 4*x3 + 12*x4  =    0 

## Example 46
A <- matrix(c(0, -9, 4, 4, 5, 6, 3, 8, 4, 6, 0, 4, 7, 7, 5, 12), 4, 4)
b <- c(-9, -17, -3, -12)
showEqn(A, b)

 0*x1 + 5*x2 + 4*x3  + 7*x4  =   -9 
-9*x1 + 6*x2 + 6*x3  + 7*x4  =  -17 
 4*x1 + 3*x2 + 0*x3  + 5*x4  =   -3 
 4*x1 + 8*x2 + 4*x3 + 12*x4  =  -12 

Working example

### Practical Applications

### Setting up library
library(matlib)

A <- matrix(c(1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 1, 0,
0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 
0, 1, 0, 0,  0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0,  0, 0,
1, 0,  0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1,
0, 0, 0, 1 ), nrow = 8, ncol = 16, byrow = TRUE)
b <- c(20, 22, 33, 21, 4, 13, 43, 36)
#b <- c(21.5, 20.9, 33.4, 22.1, 4.5, 14.2, 42.4, 62.9)
Solve(A, b)

x1         - 1*x6 - 1*x7 - 1*x8   - 1*x10 - 1*x11 - 1*x12    - 1*x14 - 1*x15 - 1*x16  =  -72 
  x2         + x6                   + x10                      + x14                  =   13 
    x3              + x7                    + x11                      + x15          =   43 
      x4                   + x8                     + x12                      + x16  =   36 
        x5   + x6   + x7   + x8                                                       =   22 
                               x9   + x10   + x11   + x12                             =   33 
                                                         x13   + x14   + x15   + x16  =   21 
                                                                                   0  =    0 

Infinitely many solutions!

Exercise