**Preface**

This series is aimed at providing tools for an electrical engineer to analyze data and solve problems in design. The focus is on applying linear algebra to systems of equations or large sets of matrix data.

**Introduction**

This article will demonstrate the use of matrix algebra to solve for equilibrium in systems. This is common in network flow, economics and electrical circuits (current and voltage analysis). We will apply exact and least squares solutions.

If you are not familiar with linear algebra or need a brush up I recommend Linear Algebra and its Applications by Lay, Lay and McDonald. It provides an excellent review of theory and applications.

This also assumes you are familiar with Python or can stumble your way through it.

The data and code are available here: linear1.py and linear1.xlsx.

**Concepts***Linear Equation:* an equation that can be written in the form *b = a _{1}x_{1} + a_{2}x_{2} + ... + a_{n}x_{n}*.

*Linear System:* A set of linear equations involving the same variables.

*Solution:* a set of values for the set of variables in a linear system that make the system a true statement.

Linear systems are either consistent (having a solution) or inconsistent (having no solution). Linear systems with a solution are either unique (1 solution) or have a free variable indicating infinitely many solutions.

*Matrix Equation:* The above definition of a linear equation uses a linear combination to create a vector equation. The same equation can be transliterated to matrix notation (used by computers) in the form **Ax** = **b**.

*Invertible Matrix Theorem:* from Lay, Lay, McDonald.

Let *A* be a square *n*x*n* matrix. Then the following statements are equivalent. That is, for a given *A*, the statements are either all true or all false.

*A*is an invertible matrix.*A*is row equivalent to the*n*x*n*identity matrix.*A*has*n*pivot positions.- The equation
*A***x**=**0**has only the trivial solution. - The columns of
*A*form a linearly independent set. - The linear transformation
**x**->*A***x**is one-to-one. - The equation
*A***x**=**b**has at least one solution for each**b**in*R*.^{n} - The columns of
*A*span*R*.^{n} - The linear transformation
**x**->*A***x**maps*R*onto^{n}*R*.^{n} - There is an
*n*x*n*matrix*C*such that*CA*=*I*. - There is an
*n*x*n*matrix*D*such that*AD*=*I*. *A*is an invertible matrix.^{T}- The columns of
*A*form a basis of*R*.^{n} - Col
*A*=*R*.^{n} - dim Col
*A*=*n* - rank
*A*=*n* - Nul
*A*={**0**} - dim Nul
*A*=0 - The number 0 is not an eigenvalue of A.
- The determinant of A is not zero.

**Importing Your Data Set**

I will use the Anaconda package suite with Python, numpy and the Spyder IDE for mathematical analysis. It is cross platform, free and open source. It is also easy to import data from Excel once you have the code snippet.

I usually have the following boilerplate code in my scripts:

import sys

import xlrd

import numpy as np

import matplotlib.pyplot as plt

from scipy.interpolate import interp1d

np.set_printoptions(threshold=np.nan)

print sys.version

print __name__

The system in this example is traffic flow for 4 intersections in a city grid. All ingress/egres traffic is known and is pictured below.

The data above has been arranged into matrix notation and stored in the excel file with the coefficients in the first five columns and the dependent variables in the sixth (last) column.

filename='linear_1.xlsx'

print

print 'Opening',filename

ss=xlrd.open_workbook(filename,on_demand=True)

sheet_index=0

print ' Worksheet',ss.sheet_names()[sheet_index]

sh=ss.sheet_by_index(sheet_index)

print ' Reading data values'

spec=[]

for row_index in range(0,sh.nrows):

spec.append([sh.cell_value(row_index,0),sh.cell_value(row_index,1),sh.cell_value(row_index,2),sh.cell_value(row_index,3),sh.cell_value(row_index,4)])

spec=np.array(spec)

print ' Design spec array size:',len(spec)

b=[]

for row_index in range(0,sh.nrows):

b.append(sh.cell_value(row_index,5))

b=np.array(b)

print ' Design b array size: ',len(b)

Set filename to your Excel filename and set sheet_index to the tab you're reading from (counting starts from zero). Set range(x,y) and {a,b,c,...,n} to the array you are reading in from the sheet. `range(0,sh.nrows)`

is a handy expression to use here.

The last bit of code converts the array to a numpy array which will be necessary for the code below.

The whole operation has fairly good performance.

**Solving the System for Equilibrium**

Now that the test results and input variables have been imported let's solve the system. We want to find the equilibrium price.

I have found that np.linalg.solve() is not always the best solution because it cannot handle systems that are underdetermined, overdetermined or linearly dependent. A better alternative seems to be the np.linalg.lstsq() function.

Both functions take the form *function(a,b)* where *a* and *b* are from *ax=b*.

#eq=np.linalg.solve(spec,b)

eq=np.linalg.lstsq(spec,b)

print

print eq

This is the output:

Opening linear_1.xlsx

Worksheet Sheet1

Reading data values

Design spec array size: 6

Design b array size: 6

(array([ 375., 425., 400., 275., 225.]), array([], dtype=float64), 4, array([ 2.11688316e+00, 1.63978419e+00, 1.41421356e+00,

9.10995891e-01, 8.65181300e-17]))

The first array is the compute result. DO NOT take that result for granted. You must examine the second array for anything unusual. The most relevant thing to note is that the fifth value is close to zero. This is a strong indication of a free (dependent) variable.

Adding a sixth row to the data set that hard sets the fifth variable to 9 results in this output:

Opening linear_1.xlsx

Worksheet Sheet1

Reading data values

Design spec array size: 6

Design b array size: 6

(array([ 591., 209., 400., 491., 9.]), array([ 4.32956906e-27]), 5, array([ 2.16121248, 1.7660084 , 1.47046231, 0.92978281, 0.42850842]))

An entirely different and this time correct result. We can see the fifth variable hard set to 9. If you plug the **x** array into each row of coefficients you'll find 100% agreement.

**Systems of Systems**

It is often desirable to solve more than one system. This can be done by iterating through each array, one at a time. This is okay but I have found that memory requirements can be utterly huge for this method. The method below seems bulky but is exponentially faster and uses less memory (I don't know how much because the last time I tried the iterative method I kept running out of memory after 4-8 hours).

The following pseudo code works well when the set of coefficients is the same and the dependent values are different. It works just as well for different sets of coefficients you just have to build your array appropriately.

I also include a determinant check to verify the array is invertible and therefore solvable.

a=np.array([[some],[set],[of],[coefficients]])

b=np.array([[some],[set],[of],[dependents]])

n=number of coefficients

if np.linalg.det(a) != 0:

csi=np.linalg.solve([a]*n,b)

**Next Up**

In the next article we will examine trend analysis using polynomial least squares.

## Add new comment