#!/usr/bin/env python
# coding: utf-8

# <style>
# @import url(https://www.numfys.net/static/css/nbstyle.css);
# </style>
# <a href="https://www.numfys.net"><img class="logo" /></a>
# 
# # Linear Algebra in Python
# 
# ### Modules - Basics
# <section class="post-meta">
# By Henning G. Hugdal, Håkon W. Ånes and Jon Andreas Støvneng
# </section>
# Last edited: February 7th 2018 
# ___

# Even though matrix operations might not be as naturally built in to Python as in for instance MATLAB, Python can be used to solve linear algebra problems. In addition to functions in the `numpy` base library, `numpy` also has its own linear algebra package called `linalg`. In this module we will take a look at how `numpy` and `linalg` can be used to handle matrices. First of all we need to import `numpy` and `linalg`:

# In[1]:


import numpy as np
import numpy.linalg as la


# ### Basic operations

# Let's first define two matrices $A$ and $B$ by using `numpy` arrays:

# In[2]:


A = np.array([[1, 2],
              [1, 2]])
B = np.array([[2, 1],
             [1, 2]])


# In order to do a regular matrix multiplication between $A$ and $B$ we use the `numpy` function `inner(A,B)`:

# In[5]:


print(np.inner(A, B))


# Powers of the same matrix can be calculated using `linalg.matrix_powers(A,n)`, where $n$ is the integer power.

# In[7]:


print(la.matrix_power(A, 2))


# ### Solving Linear Systems of Equations

# Linear systems of equations $A\vec{x}=b$ can be solved using `linalg.solve(A,b)`, which returns the vector $\vec{x}$. As an example, let's solve the system of equations
# $$\begin{align}
# 2x + y - z &= 8\\
# -3x - y + 2z &= -11\\
# -2x + y + 2z &= -3,
# \end{align}$$

# In[8]:


A = np.array([[2, 1, -1],
             [-3, -1, 2],
             [-2, 1, 2]])
b = np.array([8, -11, -3])

x = la.solve(A, b)
print(x)


#  We see that the solution is $x=2$, $y=3$ and $z=-1$. We could also have solved the above by calculating the inverse of $A$ and multiplying it with $b$. The inverse of $A$ can be calculated using `linalg.inv(A)`.

# In[9]:


A_inv = la.inv(A)
x = np.inner(A_inv, b)
print(x)


# We see that this gives the same solution, as it should!

# ### Eigenvalue Problems

# The `linalg` package can also be used to find the eigenvalus and eigenvectors of a matrix, that is, the vector $\vec{x}$ of a matrix $A$ which fulfills
# $$ A\vec{x} = a \vec{x},$$
# where $a$ is a constant. Let's determine the eigenvalues and eigenvectors of the matrix $B$ defined above.

# In[12]:


(a, x) = la.eig(B)
print("Eigenvalues:\n", a)
print("Eigenvectors:")
for i in range(len(x)):
    print(x[:, i])


# From the above we see that $B$ has the eigenvalues $3$ and $1$ with eigenvectors $[1, 1]/\sqrt{2}$ and $[-1,1]/\sqrt{2}$ respectively. The function `eig` returns an array $a$ of eigenvalues and array $x$ with normalized eigenvectors, such that the vector `x[:,i]` corresponds to the eigenvalue `a[i]`.

# The above shows just some of the various functions available in the `numpy` and `linalg` libraries. For more information, see the link to the documentation below.

# ### References
# [1] [Documentation](http://docs.scipy.org/doc/numpy/reference/routines.linalg.html) of the Numpy Linear Algebra package.