diff --git a/Tutorials/02_Linear_Algebra/2a_Linear-Algebra-Basics.ipynb b/Tutorials/02_Linear_Algebra/2a_Linear-Algebra-Basics.ipynb
new file mode 100644
index 00000000..165960a7
--- /dev/null
+++ b/Tutorials/02_Linear_Algebra/2a_Linear-Algebra-Basics.ipynb
@@ -0,0 +1,1772 @@
+{
+ "cells": [
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "\"\"\"Detailed introduction to linear algebra and matrix mechanics.\"\"\"\n",
+ "\n",
+ "__authors__ = \"D. A. Sirianni\"\n",
+ "__credits__ = [\"Daniel G. A. Smith\", \"Ashley Ringer McDonald\"]\n",
+ "__email__ = [\"sirianni.dom@gmail.com\"]\n",
+ "\n",
+ "__copyright__ = \"(c) 2014-2021, The Psi4NumPy Developers\"\n",
+ "__license__ = \"BSD-3-Clause\"\n",
+ "__date__ = \"2021-04-09\""
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Motivation & Background\n",
+ "\n",
+ "While it is possible to write down the mathematical equations which govern the physical behavior for everything\n",
+ "from electron dynamics to aeronautics, the solution of these equations is, in many cases, either too challenging\n",
+ "or actually impossible to solve exactly. It is therefore necessary to rely on methods for generating approximate\n",
+ "solutions, most of which leverage the power of computers to do so in a robust and efficient manner. While\n",
+ "computers are flexible in their ability to solve problems, they do so by representing data as discrete values --\n",
+ "either 1 or 0, true or false, on or off -- which makes solving problems from our continuous world a challenge.\n",
+ "Using a field of mathematics known as ***Linear Algebra***, however, we may transform these physcial equations\n",
+ "from their original, continuous form (most often a differential equation) to one which is amenable to being\n",
+ "solved efficiently on a computer. In doing so, the relationships between continuous variables are reframed into\n",
+ "the relationships between these continuous variables and a fixed set of discrete reference objects, referred to\n",
+ "as a _basis set_, which in turn can relate with other continuous variables. While this may sound confusing,\n",
+ "this process allows for one of the most significant advantages of computing to be applied to solving real-world\n",
+ "problems: computers can perform linear algebra operations and solve linear algebra expressions extremely\n",
+ "efficiently! \n",
+ "\n",
+ "In this lesson, we will introduce the basic principles, objects, and operations employed in linear algebra,\n",
+ "through the lens of scientific computing in the high-level Python programming language. In this way, we can\n",
+ "leverage the syntactic ease and code simplicity of Python programming to explore and experiment with these\n",
+ "linear algebra concepts, before later using them to apply quantum mechanics to describe the electronic structure\n",
+ "of atoms and molecules."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## The Basis for Decretizing Continuous Variables\n",
+ "\n",
+ "In life, we often note the movement of objects --- even ourselves --- in relation to other objects, each of\n",
+ "which may or may not be moving too. For instance, in order to describe the manner in which every object in the\n",
+ "universe moves relative to every other object in the universe is a near-infinitely complicated problem. To\n",
+ "simplify the situation somewhat, let's ignore the \"rest of the universe\" other than your immediate vicinity,\n",
+ "e.g., the room you currently occupy. Every object in the room, yourself included, exerts a gravitational pull on\n",
+ "every other object in the room, according to Newton's universal law of gravitation. To fully describe these \n",
+ "gravitational forces, therefore, it would be seemingly necessary to keep track not only of every object's position\n",
+ "and movement relative to the other objects, but also _every object's position and movement relative to every other\n",
+ "object._ Even in your immediate vicinity, that is a lot of information! Shouldn't there be some easier way to\n",
+ "keep track of position or movement?\n",
+ "\n",
+ "Let's start from our own perspective. Some short distance away from where you are, is the screen on which you\n",
+ "are reading these words. More than with just this linear distance, however, we can describe the position of the\n",
+ "screen relative to your eyes in terms of its _vertical displacement_ (i.e., how far up/down you are looking),\n",
+ "its _lateral displacement_ (i.e., how close/far into the distance you must focus your eyes), and its _horizontal\n",
+ "displacement_ (i.e., how far left/right on the screen your eyes are). Breaking down the screen position relative\n",
+ "to your eyes into these up/down, close/far, and left/right _components_ is precisely the manner in which we can\n",
+ "simplify the definition of our surroundings: rather than defining every object's position relative to each other,\n",
+ "we may do so by defining each object's position relative to a single fixed object called the _origin_ and fixed\n",
+ "directions which form the _basis_ for our definitions of position from the origin. While defining a particular\n",
+ "origin and basis should depend on what is most sensible in a given scenario, these are sufficiently general\n",
+ "concepts that we may then use to develop a practical framework for linear algebra.\n",
+ "\n",
+ "### What is a vector?\n",
+ "\n",
+ "Unlike ordinary numbers, e.g., 1, -23.442, $\\sqrt{5}$, etc. which have only a magnitude (and which we will refer\n",
+ "to as _scalars_), vectors are mathematical quantities with both a _magnitude_ and a _direction_. For example,\n",
+ "the distance someone walks (say, 3 miles) is a scalar quantity, but we could make this into a vector by adding\n",
+ "the information that the person walks 3 miles due North. Denoting vectors in terms of a standard set of reference\n",
+ "directions is common practice -- in fact, we've just used the cardinal directions (North, South, East, West) to\n",
+ "do so. But what about other kinds of vectors? What about when someone throws a ball at a 35$^{\\circ}$ angle above\n",
+ "the horizontal? How would we describe the direction of that vector?\n",
+ "\n",
+ "Most often, vectors are denoted as an ordered collection of scalar _components_ which describe the magnitude\n",
+ "of the vector in the directions of several _basis vectors_. In our example above, if we define the North-South\n",
+ "line to be the positive and negative $y$-axis, and similarly for East-West to be the positive and negative\n",
+ "$x$-axis, then the vector describing a person walking 3 miles due north could be represented as an ordered pair of\n",
+ "$x$ and $y$ components:\n",
+ "\n",
+ "$${\\bf v} = \\begin{pmatrix} 0 & 3\\end{pmatrix} = \\begin{pmatrix} v_x & v_y\\end{pmatrix}.$$\n",
+ "\n",
+ "Here, the vectors ${\\bf e_1} = \\begin{pmatrix} 1 & 0\\end{pmatrix}$ and ${\\bf e_2} = \\begin{pmatrix} 0 & 1\n",
+ "\\end{pmatrix}$, each running along the $x$ and $y$ axes, respectively, forms the _basis_ within which we define\n",
+ "the vector ${\\bf v}$, and $v_x = 0$, $v_y = 3$ are the components of ${\\bf v}$ in each of these basis vectors'\n",
+ "directions. While this example used a vector which has two components, $v_x$ and $v_y$, vectors can have any\n",
+ "number of components, and is more generally represented as\n",
+ "\n",
+ "$${\\bf v} = \\begin{pmatrix} v_1 & v_2 & \\cdots & v_n\\end{pmatrix},$$\n",
+ "\n",
+ "which we say has _length_ $n$ because it has $n$ components. In fact, the movement of a baseball as it is thrown\n",
+ "would be best described by a length-3 vector, with components in each of the $x$, $y$, and $z$ directions, and\n",
+ "we will see in the future that in computational chemistry, vectors can easily have lengths in the millions or even\n",
+ "_billions_. Bet you're glad we're using a computer to do that work, huh?\n",
+ "\n",
+ "### Representing Vectors in Python\n",
+ "\n",
+ "So far, we have learned that a vector is simply an ordered collection of scalar components. Therefore, we can\n",
+ "represent these objects with any Python type that is a similarly ordered collection, including a `list` or \n",
+ "`tuple`. In the cell below, we will first define two length-3 vectors ${\\bf v}$ and ${\\bf w}$:\n",
+ "\n",
+ "\\begin{align}\n",
+ "{\\bf v} = \\begin{pmatrix} 1 & 2 & 3\\end{pmatrix}\\\\\n",
+ "{\\bf w} = \\begin{pmatrix} 4 & 5 & 6\\end{pmatrix}\\\\\n",
+ "\\end{align}"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# ==> Representing Vectors <==\n",
+ "# Define two length-3 vectors, v & w\n",
+ "v = (1, 2, 3) # Define as a tuple\n",
+ "w = [4, 5, 6] # Define as a list"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "While both `list` and `tuple` types seem perfectly adequate to represent the ordered structure of vectors,\n",
+ "they behave very differently in practice and generally should not be mixed. To illustrate this difference\n",
+ "and to see how this could be problematic, let's say that we made a mistake when we defined our vectors above,\n",
+ "where each element should actually be scaled by a factor of 10, i.e., they really should be defined as\n",
+ "\n",
+ "\\begin{align}\n",
+ "{\\bf v} = \\begin{pmatrix} 10 & 20 & 30\\end{pmatrix}\\\\\n",
+ "{\\bf w} = \\begin{pmatrix} 40 & 50 & 60\\end{pmatrix}.\n",
+ "\\end{align}\n",
+ "\n",
+ "Execute the cells below to update the values of each element of our two vectors to be 10 times larger,\n",
+ "using `for` loops."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 3,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[40, 50, 60]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Redefine elements of `w` using a `for` loop <==\n",
+ "for i in range(len(w)):\n",
+ " w[i] *= 10\n",
+ " \n",
+ "print(w)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 4,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Tuples are _immutable_, so you can't change their elements after creation!\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Try to redefine elements of `v` using a `for` loop <==\n",
+ "try:\n",
+ " for i in range(len(v)):\n",
+ " v[i] *= 10\n",
+ "\n",
+ " print(v)\n",
+ "\n",
+ "except TypeError:\n",
+ " print(\"Tuples are _immutable_, so you can't change their elements after creation!\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Uh-oh! \n",
+ "\n",
+ "The reassignment of the elements of `w` seemed to work just fine, but the same approach for `v` failed\n",
+ "with a `TypeError`. This is because unlike `list`s, `tuple`s are _immutable_ types: in other words, once a vector\n",
+ "is created as a tuple, it can never be changed in any way. Unfortunately for us, that would mean that\n",
+ "the rest of the lesson -- where we finally get to _do_ fun things with our vectors and matrices -- would be\n",
+ "pointless, because our objects could never change! Therefore until we begin to use specialized data types\n",
+ "specifically designed to represent arrays, we will stick with `list`s to define our vectors and matrices.\n",
+ "\n",
+ "To correct this problem so that we may continue the lesson without encountering `tuple`-related `TypeError`s,\n",
+ "use the cell below to first redefine the vector `v` as a list, before then updating its values such that\n",
+ "\n",
+ "$${\\bf v} = \\begin{pmatrix} 10 & 20 & 30\\end{pmatrix}.$$"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 5,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[10, 20, 30]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Redefine v <==\n",
+ "# Redefine as list\n",
+ "v = [1, 2, 3] # Could do directly\n",
+ "v = list(v) # Could also do w/ typecasting\n",
+ "\n",
+ "# Update elements of v using for-loop\n",
+ "for i in range(len(v)):\n",
+ " v[i] *= 10\n",
+ " \n",
+ "print(v)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Vector Operations\n",
+ "\n",
+ "Now that we know what vectors are and how to properly represent them in Python, we can begin to actually _do_ something\n",
+ "with them other than just storing values inside of them! As it turns out, vectors can interact with scalar values\n",
+ "and each other through some of the same operations that scalar values interact with each other, namely addition\n",
+ "and multiplication. Because of the additional structure that vectors possess, however, these operations are\n",
+ "slightly more complicated than for pure scalars. To introduce these new vector operations, we will use the\n",
+ "vectors **v** and **w** we created above.\n",
+ "\n",
+ "### Vector Addition\n",
+ "For two vectors ${\\bf v} = \\begin{pmatrix} v_1 & v_2 & \\cdots & v_n \\end{pmatrix}$ and\n",
+ "${\\bf w} = \\begin{pmatrix} w_1 & w_2 & \\cdots & w_n \\end{pmatrix}$,\n",
+ "\n",
+ "\n",
+ "$${\\bf z} = {\\bf v} + {\\bf w} = \\begin{pmatrix} v_1 & v_2 & \\cdots & v_n \\end{pmatrix} +\n",
+ "\\begin{pmatrix} w_1 & w_2 & \\cdots & w_n \\end{pmatrix} = \n",
+ "\\begin{pmatrix} v_1 + w_1 & v_n + w_2 & \\cdots & v_n + w_n \\end{pmatrix}\n",
+ "$$\n",
+ "\n",
+ "In the cell below, evaluate the vector sum ${\\bf z} = {\\bf v} + {\\bf w}$, using the vectors we defined above\n",
+ "and a `for`-loop:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 6,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[50, 70, 90]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Define function to evaluate vector sum v + w using a for loop, storing result in z <==\n",
+ "\n",
+ "def vector_add(v, w):\n",
+ " \n",
+ " z = [0] * len(v)\n",
+ " \n",
+ " for i in range(len(v)):\n",
+ " z[i] = v[i] + w[i]\n",
+ " \n",
+ " return z\n",
+ "\n",
+ "print(vector_add(v,w))"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Notice that **z**, the sum of **v** + **w**, is the same length as both **v** and **w** themselves. Consequently,\n",
+ "vector addition requires that the vectors being added have the same length.\n",
+ "\n",
+ "### Scalar Addition & Multiplication for Vectors\n",
+ "\n",
+ "For a vector ${\\bf v} = \\begin{pmatrix} v_1 & v_2 & \\cdots & v_n \\end{pmatrix}$ and\n",
+ "scalar (i.e., regular numbers) values $r$ and $s$, we define _scalar multiplication_ and _scalar addition_ as\n",
+ "\n",
+ "$${\\bf z} = r \\cdot {\\bf v} + s = r \\cdot \\begin{pmatrix} v_1 & v_2 & \\cdots & v_n \\end{pmatrix} + s = \n",
+ "\\begin{pmatrix} r \\cdot v_1 + s & r \\cdot v_2 + s & \\cdots & r \\cdot v_n + s \\end{pmatrix}\n",
+ "$$\n",
+ "\n",
+ "In the cell below, evaluate the expression $z = r \\cdot {\\bf v} + s$, using the vector **v** we defined above\n",
+ "and the scalars r=3, s = 1:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 7,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[31, 61, 91]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Scalar Multiplication & Addition <==\n",
+ "def rvps(v, r, s):\n",
+ " z = [0] * len(v)\n",
+ " \n",
+ " for i in range(len(v)):\n",
+ " z[i] = r * v[i] + s\n",
+ " \n",
+ " return z\n",
+ "\n",
+ "print(rvps(v, 3, 1))"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Vector Multiplication\n",
+ "\n",
+ "Unlike with scalars, there exist several ways to perform vector multiplication, due to the added structure\n",
+ "of vectors. \n",
+ "\n",
+ "#### Elementwise Vector Product\n",
+ "\n",
+ "The most straightforward product for two vectors **v** and **w** is the _elementwise_ product,\n",
+ "which we denote with the $\\odot$ symbol (`*` when in code), given by:\n",
+ "\n",
+ "$$ {\\bf z} = {\\bf v}\\odot{\\bf w} = \\begin{pmatrix} v_1 & v_2 & \\cdots & v_n \\end{pmatrix} \\odot\n",
+ "\\begin{pmatrix} w_1 & w_2 & \\cdots & w_n \\end{pmatrix} = \n",
+ "\\begin{pmatrix} v_1 \\cdot w_1 & v_n \\cdot w_2 & \\cdots & v_n \\cdot w_n \\end{pmatrix}\n",
+ "$$\n",
+ "\n",
+ "Using **v** and **w** defined above, compute their elementwise product, ${\\bf v}\\odot{\\bf w}$, in the cell below."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 8,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[400, 1000, 1800]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Elementwise product, v * w <==\n",
+ "\n",
+ "def vector_lmntproduct(v, w):\n",
+ " \n",
+ " z = [0] * len(v)\n",
+ " \n",
+ " for i in range(len(v)):\n",
+ " z[i] = v[i] * w[i]\n",
+ " \n",
+ " return z\n",
+ "\n",
+ "print(vector_lmntproduct(v,w))"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "While the elementwise vector product is used in image compression and machine learning, it is less\n",
+ "useful in physics-based applications than other types of vector multiplication, which we will explore\n",
+ "below.\n",
+ "\n",
+ "#### Vector Dot Product\n",
+ "\n",
+ "For our vectors ${\\bf v}$ and ${\\bf w}$, the dot product ${\\bf z} = {\\bf v}\\cdot{\\bf w}$ is given by\n",
+ "\n",
+ "$$\n",
+ "{\\bf z} = {\\bf v}\\cdot{\\bf w} = \\begin{pmatrix} v_1 & v_2 & \\cdots & v_n \\end{pmatrix} \\cdot\n",
+ "\\begin{pmatrix} w_1 & w_2 & \\cdots & w_n \\end{pmatrix} = v_1\\cdot w_1 + v_2\\cdot w_2 + \\ldots + v_n\\cdot w_n\n",
+ "$$\n",
+ "\n",
+ "Notice that the dot product of two vectors actually yields a scalar, rather than another vector. This scalar\n",
+ "value has special importance relating ${\\bf v}$ and ${\\bf w}$, since\n",
+ "\n",
+ "$${\\bf z} = {\\bf v}\\cdot{\\bf w} = \\vert{\\bf v}\\vert\\cdot\\vert{\\bf w}\\vert\\cos{\\theta},$$\n",
+ "\n",
+ "where $\\theta$ is the angle between the vectors ${\\bf v}$ and ${\\bf w}$. Therefore, the dot product is a measure\n",
+ "of the _overlap_ between two vectors, or the extent to which the vectors have the same direction and magnitude.\n",
+ "\n",
+ "In the cell below, write a function to evaluate the dot product between two vectors, and use it to evaluate\n",
+ "${\\bf v}\\cdot{\\bf w}$, ${\\bf w}\\cdot {\\bf v}$, and ${\\bf v}\\cdot {\\bf v}$.\n",
+ "\n",
+ "> Note: We denote the dot product ${\\bf v}\\cdot {\\bf w}$ as `< v | w >` in code comments to differentiate it\n",
+ "from the elementwise product, ${\\bf v}\\odot{\\bf w}$ (`v * w` in code)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 9,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "< v | w > = [400, 1000, 1800]\n",
+ "< w | v > = [400, 1000, 1800]\n",
+ "< v | v > = [100, 400, 900]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Dot product practice <==\n",
+ "# Define general dot product function\n",
+ "def vector_dot(v, w):\n",
+ " z = list(range(len(v)))\n",
+ " # Check lengths of v & w are equal\n",
+ " assert len(v)==len(w), f\"Vector arguments do not have equal length!\"\n",
+ " # Compute dot product\n",
+ " for i in range(len(v)):\n",
+ " z[i] = v[i] * w[i]\n",
+ " \n",
+ " return z\n",
+ "\n",
+ "# Evaluate < v | w >\n",
+ "print(f\"< v | w > = {vector_dot(v, w)}\")\n",
+ "\n",
+ "# Evaluate < w | v >\n",
+ "print(f\"< w | v > = {vector_dot(w, v)}\")\n",
+ "\n",
+ "# Evaluate v^2 = < v | v >\n",
+ "print(f\"< v | v > = {vector_dot(v, v)}\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "#### Vector Cross Product\n",
+ "\n",
+ "Finally, the _cross product_ of vectors ${\\bf v}$ and ${\\bf w}$, denoted ${\\bf v}\\times{\\bf w}$, is given by\n",
+ "\n",
+ "$${\\bf z} = {\\bf v}\\times{\\bf w} = \\begin{pmatrix} v_2w_3 - w_2v_3 & v_1w_3 - w_1v_3 & v_1w_2 - w_1v_2\\end{pmatrix},$$\n",
+ "\n",
+ "which is a vector perpendicular to both ${\\bf v}$ and ${\\bf w}$. While the cross product of vectors is\n",
+ "exceptionally useful in classical physical theories, particularly in Maxwell's formulation of electromagnetism,\n",
+ "we will not generally use the cross product in computational chemistry applications. "
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## What is a Matrix?\n",
+ "\n",
+ "Just like a vector is an ordered 1-dimensional collection of scalars (i.e., the scalars occupy a single row), \n",
+ "a _matrix_ is a 2-dimensional ordered collection of scalars, organized into rows _and_ columns:\n",
+ "\n",
+ "$$\n",
+ "{\\bf M} = \\begin{pmatrix}\n",
+ "M_{11} & M_{12} & \\cdots & M_{1n}\\\\\n",
+ "M_{21} & M_{22} & \\cdots & M_{2n}\\\\\n",
+ "\\vdots & \\vdots & \\ddots & \\vdots\\\\\n",
+ "M_{m1} & M_{m2} & \\cdots & M_{mn}\n",
+ "\\end{pmatrix}\n",
+ "$$\n",
+ "\n",
+ "Here, ${\\bf M}$ is a $m\\times n$ matrix, and in general $m$ does not have to equal $n$, i.e., ${\\bf M}$ does\n",
+ "not have to be _square_. One useful way to think of matrices is as an ordered collection of _vectors_:\n",
+ "\n",
+ "$$\n",
+ "{\\bf M} = \\begin{pmatrix}\n",
+ "M_{11} & M_{12} & \\cdots & M_{1n}\\\\\n",
+ "M_{21} & M_{22} & \\cdots & M_{2n}\\\\\n",
+ "\\vdots & \\vdots & \\ddots & \\vdots\\\\\n",
+ "M_{m1} & M_{m2} & \\cdots & M_{mn}\n",
+ "\\end{pmatrix} = \\begin{pmatrix}\n",
+ "{\\bf v}_1\\\\\n",
+ "{\\bf v}_2\\\\\n",
+ "\\vdots\\\\\n",
+ "{\\bf v}_m\n",
+ "\\end{pmatrix},\n",
+ "$$\n",
+ "\n",
+ "where ${\\bf v}_i = \\begin{pmatrix}M_{i1} & M_{i2} & \\cdots & M_{in}\\end{pmatrix}$ is the $i$th _row vector_\n",
+ "of ${\\bf M}$. \n",
+ "\n",
+ "### Representing Matrices in Python\n",
+ "\n",
+ "Since a matrix can be thought of as an ordered collection of its rows, it seems sensible to represent a matrix\n",
+ "as a `list` of `list`s. Using this principle, in the cell below define the following matrix:\n",
+ "\n",
+ "$${\\bf A} = \\begin{pmatrix}\n",
+ "1 & 2 & 3\\\\\n",
+ "4 & 5 & 6\\\\\n",
+ "7 & 8 & 9\n",
+ "\\end{pmatrix}\n",
+ "$$"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 10,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[[1, 2, 3], [4, 5, 6], [7, 8, 9]]\n",
+ "[[10, 11, 12], [13, 14, 15], [16, 17, 18]]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Define matrices A and B as list of lists <==\n",
+ "\n",
+ "A = [[1, 2, 3], [4, 5, 6], [7, 8, 9]]\n",
+ "B = [[10, 11, 12], [13, 14, 15], [16, 17, 18]]\n",
+ "\n",
+ "print(A)\n",
+ "print(B)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "While it is certainly possible to represent a matrix as a collection of rows, it is equally possible to define\n",
+ "a matrix as a collection of columns:\n",
+ "\n",
+ "$$\n",
+ "{\\bf M} = \\begin{pmatrix}\n",
+ "M_{11} & M_{12} & \\cdots & M_{1n}\\\\\n",
+ "M_{21} & M_{22} & \\cdots & M_{2n}\\\\\n",
+ "\\vdots & \\vdots & \\ddots & \\vdots\\\\\n",
+ "M_{m1} & M_{m2} & \\cdots & M_{mn}\n",
+ "\\end{pmatrix} = \\begin{pmatrix}\n",
+ "{\\bf v}_1\\\\\n",
+ "{\\bf v}_2\\\\\n",
+ "\\vdots\\\\\n",
+ "{\\bf v}_m\n",
+ "\\end{pmatrix} = \\begin{pmatrix}\n",
+ "{\\bf w}_1 & {\\bf w}_2 & \\cdots & {\\bf w}_n\n",
+ "\\end{pmatrix},\n",
+ "$$\n",
+ "\n",
+ "where as above, ${\\bf v}_i = \\begin{pmatrix}M_{i1} & M_{i2} & \\cdots & M_{in}\\end{pmatrix}$ is the $i$th row\n",
+ "vector, but now ${\\bf w}_j = \\begin{pmatrix}M_{1j} \\\\ M_{2j} \\\\ \\vdots \\\\ M_{mj}\\end{pmatrix}$ is the $j$th\n",
+ "_column vector_ of ${\\bf M}$. Now, representing a matrix as a `list` of `list`s seems less straightforward,\n",
+ "since in Python there is no difference between row and column vectors (they're both just `list`s!).\n",
+ "So, when constructing a matrix from a collection of vectors represented as `list`s, the final matrix will\n",
+ "depend upon whether each vector is assumed to represent a row or column of the matrix. To illustrate this,\n",
+ "consider the $3\\times 3$ matrix \n",
+ "\n",
+ "$${\\bf M} = \\begin{pmatrix}\n",
+ "1 & 2 & 3\\\\\n",
+ "1 & 2 & 3\\\\\n",
+ "1 & 2 & 3\n",
+ "\\end{pmatrix}\n",
+ "$$\n",
+ "\n",
+ "By defining the row vectors `r_i = [1, 2, 3]` and column vectors `c_i = [i, i, i]`, try to form the matrix\n",
+ "`M` which matches the one above."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 11,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[[1, 2, 3], [1, 2, 3], [1, 2, 3]]\n",
+ "[[1, 1, 1], [2, 2, 2], [3, 3, 3]]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Matrix formation from row vectors vs. column vectors <==\n",
+ "\n",
+ "# Define row (r_i) & column vectors (c_i)\n",
+ "r_1 = [1, 2, 3]\n",
+ "r_2 = [1, 2, 3]\n",
+ "r_3 = [1, 2, 3]\n",
+ "\n",
+ "Mrow = [r_1,\n",
+ " r_2,\n",
+ " r_3]\n",
+ "\n",
+ "# Try to form M as a vector of columns or rows\n",
+ "c_1 = [1, 1, 1]\n",
+ "c_2 = [2, 2, 2]\n",
+ "c_3 = [3, 3, 3]\n",
+ "\n",
+ "Mcol = [c_1, c_2, c_3]\n",
+ "\n",
+ "print(Mrow)\n",
+ "print(Mcol)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 12,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[[1, 2, 3], [1, 2, 3], [1, 2, 3]]\n",
+ "[[1, 1, 1], [2, 2, 2], [3, 3, 3]]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Matrix formation from row/column vectors <==\n",
+ "\n",
+ "# Define row (r_i) & column vectors (c_i)\n",
+ "r_1 = [1, 2, 3]\n",
+ "r_2 = [1, 2, 3]\n",
+ "r_3 = [1, 2, 3]\n",
+ "Mrow = [r_1,\n",
+ " r_2,\n",
+ " r_3]\n",
+ "\n",
+ "# Try to form M as a vector of columns or rows\n",
+ "c_1 = [1, 1, 1]\n",
+ "c_2 = [2, 2, 2]\n",
+ "c_3 = [3, 3, 3]\n",
+ "Mcol = [c_1, c_2, c_3]\n",
+ "\n",
+ "print(Mrow)\n",
+ "print(Mcol)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "As you can see, the `Mrow` matrix reproduces our definition of ${\\bf M}$:\n",
+ "\n",
+ "$${\\bf M}_{\\rm row} = \\begin{pmatrix} {\\bf r}_1\\\\ {\\bf r}_2\\\\ {\\bf r}_3\\end{pmatrix} = \n",
+ "\\begin{pmatrix}\n",
+ "1 & 2 & 3\\\\\n",
+ "1 & 2 & 3\\\\\n",
+ "1 & 2 & 3\n",
+ "\\end{pmatrix} = {\\bf M},\n",
+ "$$\n",
+ "\n",
+ "while `Mcol` does not:\n",
+ "\n",
+ "$${\\bf M}_{\\rm col} = \\begin{pmatrix} {\\bf c}_1 & {\\bf c}_2 & {\\bf c}_3\\end{pmatrix} = \n",
+ "\\begin{pmatrix}\n",
+ "1 & 1 & 1\\\\\n",
+ "2 & 2 & 2\\\\\n",
+ "3 & 3 & 3\n",
+ "\\end{pmatrix} \\neq {\\bf M}.\n",
+ "$$\n",
+ "\n",
+ "Clearly, by using `list`s to represent matrices in Python, we have implicitly assumed that matrices are formed\n",
+ "by row vectors, since `Mrow` matches ${\\bf M}$ above, but `Mcol` does not. Even though `Mcol` is not identical\n",
+ "to ${\\bf M}$ and `Mrow`, however, the two matrices do seem to be related somehow...\n",
+ "\n",
+ "### Matrix Operations\n",
+ "\n",
+ "#### Matrix Transpose\n",
+ "\n",
+ "The relationship between the matrices `Mcol` and `Mrow` is referred to as the _matrix transpose_, which is\n",
+ "a _unary_ matrix operation. In contrast to a _binary_ operation (like addition or multiplication) \n",
+ "which modifies two objects, a unary operation modifies only a single object. For a general $M\\times N$ matrix\n",
+ "${\\bf M}$, the transpose of ${\\bf M}$, ${\\bf M}^{\\rm T}$, is obtained by switching its rows and columns:\n",
+ "\n",
+ "$$\n",
+ "{\\bf M}^{\\rm T} = \\begin{pmatrix}\n",
+ "M_{11} & M_{12} & \\cdots & M_{1n}\\\\\n",
+ "M_{21} & M_{22} & \\cdots & M_{2n}\\\\\n",
+ "\\vdots & \\vdots & \\ddots & \\vdots\\\\\n",
+ "M_{m1} & M_{m2} & \\cdots & M_{mn}\n",
+ "\\end{pmatrix}^{\\rm T} = \\begin{pmatrix}\n",
+ "M_{11} & M_{21} & \\cdots & M_{m1}\\\\\n",
+ "M_{12} & M_{22} & \\cdots & M_{m2}\\\\\n",
+ "\\vdots & \\vdots & \\ddots & \\vdots\\\\\n",
+ "M_{1n} & M_{2n} & \\cdots & M_{nm}\n",
+ "\\end{pmatrix}\n",
+ "$$\n",
+ "\n",
+ "In the cell below, define a function to return the transpose of a rectangular matrix, and use it to verify\n",
+ "that `Mrow` and `Mcol` are indeed related via the transpose."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 13,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[[1, 2, 3], [1, 2, 3], [1, 2, 3]]\n",
+ "[[1, 2, 3], [1, 2, 3], [1, 2, 3]]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Define matrix transpose <==\n",
+ "def transpose(M):\n",
+ " # Get shape of M: m x n\n",
+ " m = len(M) # Number of rows\n",
+ " n = len(M[0]) # Number of columns\n",
+ " \n",
+ " # Define n x m zero matrix for M.T\n",
+ " MT = [[0 * j for j in range(m)] for i in range(n)]\n",
+ " \n",
+ " # Swap rows and columns in M to populate MT\n",
+ " for i in range(n):\n",
+ " for j in range(m):\n",
+ " MT[i][j] = M[j][i]\n",
+ " \n",
+ " return MT\n",
+ "\n",
+ "# Verify Mrow = Mcol.T\n",
+ "print(Mrow)\n",
+ "print(transpose(Mcol))"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "So, we see the reason that `Mcol` and `Mrow` were not equivalent is because matrices represented as `list`s assume\n",
+ "that the component vectors are row vectors, _which are themselves the transpose of column vectors_: \n",
+ "\n",
+ "\\begin{align}\n",
+ "{\\bf v}_{\\rm row} &= \\begin{pmatrix} v_1 & v_2 & \\cdots & v_n\\end{pmatrix} = \n",
+ "\\begin{pmatrix} v_1\\\\ v_2\\\\ \\vdots\\\\ v_n\\end{pmatrix}^{\\rm T} = {\\bf v}_{\\rm column}^{\\rm T}\\\\\n",
+ "{\\bf v}_{\\rm column} &= \\begin{pmatrix} v_1\\\\ v_2\\\\ \\vdots\\\\ v_n\\end{pmatrix} = \n",
+ "\\begin{pmatrix} v_1 & v_2 & \\cdots & v_n\\end{pmatrix}^{\\rm T} = {\\bf v}_{\\rm row}^{\\rm T}\\\\\n",
+ "\\end{align}\n",
+ "\n",
+ "#### Row Space vs. Column Space\n",
+ "\n",
+ "While it may seem like the discussion of row vectors vs. column vectors and their relationship via the transpose\n",
+ "operation was an unnecessary diversion, this turns out to be a very important concept. We can even go one step\n",
+ "further, to consider vector spaces (like our 3-dimensional world) which are defined using row vectors as being\n",
+ "distinct from those defined using column vectors; we will distinguish such vector spaces by referring to them as \n",
+ "either a _row space_, defined by row vectors, or as a _column space_ defined by column vectors. By default, we\n",
+ "will assume that all 1-dimensional arrays are column vectors, and therefore that we are working within a column\n",
+ "space.\n",
+ "\n",
+ "### Binary Matrix Operations\n",
+ "\n",
+ "#### Matrix Addition\n",
+ "For matrices **A** and **B**, we define _matrix addition_ as\n",
+ "$${\\bf C} = {\\bf A} + {\\bf B} = \\begin{pmatrix}\n",
+ "a & b\\\\\n",
+ "c & d\n",
+ "\\end{pmatrix} + \\begin{pmatrix}\n",
+ "e & f\\\\\n",
+ "g & h\n",
+ "\\end{pmatrix} = \\begin{pmatrix}\n",
+ "a + e & b + f\\\\\n",
+ "c + g & d + h\n",
+ "\\end{pmatrix}\n",
+ "$$\n",
+ "\n",
+ "In the cell below, write a function to add the matrices **A** and **B** we defined above using `for` loops:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 14,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[[11, 13, 15], [17, 19, 21], [23, 25, 27]]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Implement C = A + B using for loops <==\n",
+ "def matrix_add(A, B):\n",
+ " # Get shape of A: Ar x Ac\n",
+ " Ar = len(A) # Number of rows\n",
+ " Ac = len(A[0]) # Number of columns\n",
+ " \n",
+ " # Define Ar x Ac zero matrix to store A + B\n",
+ " C = [[0 * j for j in range(Ac)] for i in range(Ar)]\n",
+ " \n",
+ " # Compute the matrix addition & populate C with a double-for-loop\n",
+ " for i in range(len(A)):\n",
+ " for j in range(len(A[i])):\n",
+ " C[i][j] = A[i][j] + B[i][j]\n",
+ " \n",
+ " return C\n",
+ " \n",
+ "print(matrix_add(A, B))"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "#### Scalar Multiplication & Addition\n",
+ "For a matrix **A** and scalars _r_ and _s_, we define scalar multiplication and addition as\n",
+ "\n",
+ "$$r\\cdot{\\bf A} + s= r\\cdot\\begin{pmatrix}\n",
+ "a & b\\\\\n",
+ "c & d\n",
+ "\\end{pmatrix} + s = \\begin{pmatrix}\n",
+ "r\\cdot a + s & r\\cdot b + s \\\\\n",
+ "r\\cdot c + s & r\\cdot d + s \n",
+ "\\end{pmatrix}\n",
+ "$$\n",
+ "\n",
+ "In the cell below, write another function using a `for` loop to evaluate $r{\\bf A} + s$, with **A** defined\n",
+ "above and $r=2,\\,s=5$."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 15,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[[7, 9, 11], [13, 15, 17], [19, 21, 23]]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Implement C = r * A + s using for loops <==\n",
+ "def rAps(A, r, s):\n",
+ " # Get shape of A: Ar x Ac\n",
+ " Ar = len(A) # Number of rows\n",
+ " Ac = len(A[0]) # Number of columns\n",
+ " \n",
+ " # Define Ar x Ac zero matrix to store A + B\n",
+ " C = [[0 * j for j in range(Ac)] for i in range(Ar)]\n",
+ " \n",
+ " # Compute the r*A + s & populate C with a double-for-loop\n",
+ " for i in range(len(A)):\n",
+ " for j in range(len(A[i])):\n",
+ " C[i][j] = r * A[i][j] + s\n",
+ " \n",
+ " return C\n",
+ " \n",
+ "print(rAps(A, r=2, s=5))"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "#### Matrix Multiplication\n",
+ "\n",
+ "Matrix multiplication is slightly trickier than matrix addition, but has a simple pattern:\n",
+ "\n",
+ "In other words, the $i,k$-th entry of the product array is the vector dot product\n",
+ "\n",
+ "$${\\bf C} = {\\bf A}\\times{\\bf B} = \\sum_{i=1}^{M}\\sum_{k=1}^{N}\\sum_{j=1}^{P}A_{ik}B_{kj}$$\n",
+ "\n",
+ "One caveat to this matrix-matrix multiplication is that, like other matrix operations, the arrays must have\n",
+ "compatible shapes. In the case of matrix multiplication, the _inner dimensions_ of the matrices must be equal:\n",
+ "i.e., a $5\\times 2$ matrix can be multiplied by a $2\\times 4$ matrix, but not by a $3\\times 4$ matrix. If two\n",
+ "matrices are compatible, their matrix product will then have the shape of the _outer dimensions_ of the input\n",
+ "arrays, i.e., a $5\\times 2$ matrix multiplied by a $2\\times 4$ matrix will yield a $5\\times 4$ matrix.\n",
+ "\n",
+ "\n",
+ "In the cell below, define a function to return the product of two matrices, making sure to check that they\n",
+ "are compatible:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 16,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[[84, 90, 96], [201, 216, 231], [318, 342, 366]]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Implement the matrix product of A x B using Python for-loops <==\n",
+ "\n",
+ "def MM(A, B):\n",
+ " # Get shape of A: Ar x Ac\n",
+ " Ar = len(A) # Number of rows\n",
+ " Ac = len(A[0]) # Number of columns\n",
+ " \n",
+ " # Get shape of B: Br x Bc\n",
+ " Br = len(B) # Number of rows\n",
+ " Bc = len(B[0]) # Number of columns\n",
+ " \n",
+ " # Are A & B compatible? Use an assert statement to check that \"inner\" dimensions match\n",
+ " assert Ac == Br, f\"Matrices {A} and {B} are not compatible for matrix multiplication\"\n",
+ " \n",
+ " # Define Ar x Bc zero matrix to store A + B\n",
+ " C = [[0 * j for j in range(Bc)] for i in range(Ar)]\n",
+ " \n",
+ " # Evaluate AxB & populate C using a triple-for-loop\n",
+ " for i in range(len(C)):\n",
+ " for j in range(len(C[i])):\n",
+ " for k in range(len(B)):\n",
+ " C[i][j] += A[i][k] * B[k][j]\n",
+ " return C\n",
+ " \n",
+ "\n",
+ "print(MM(A, B))\n"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Redefining Vector Products as Matrix Products\n",
+ "\n",
+ "As we will see below, matrix products may be generalized to be applicable to arrays which are three-, four-, and\n",
+ "multi-dimensional. In the same way, we may actually write the vector products introduced above as matrix products.\n",
+ "As we will see, writing vector products in this manner not only conveys additional information, but also\n",
+ "illuminates new operations and provides additional flexibility.\n",
+ "\n",
+ "#### Dot Product of Column/Row Vectors\n",
+ "\n",
+ "We may redefine the simple product defined above can as:\n",
+ "\n",
+ "\\begin{align}\n",
+ "{\\bf v}\\cdot{\\bf w} = \\sum_i v_i w_i &= {\\bf v}_{\\rm row}{\\bf w}_{\\rm row}^{\\rm T} = \\begin{pmatrix} v_1 & v_2 & \\cdots & v_n\\end{pmatrix} \\cdot \\begin{pmatrix}w_1 & w_2 & \\cdots & w_n\\end{pmatrix}^{\\rm T} = \\begin{pmatrix}v_1 & v_2 & \\cdots & v_n\\end{pmatrix}\\begin{pmatrix} w_1\\\\ w_2\\\\ \\vdots\\\\ w_n\\end{pmatrix}\\\\\n",
+ "&= {\\bf v}_{\\rm col}^{\\rm T}{\\bf w}_{\\rm col} = \\begin{pmatrix} v_1 \\\\ v_2 \\\\ \\vdots \\\\ v_n\\end{pmatrix}^{\\rm T} \\cdot \\begin{pmatrix}w_1 \\\\ w_2 \\\\ \\vdots \\\\ w_n\\end{pmatrix} = \\begin{pmatrix}v_1 & v_2 & \\cdots & v_n\\end{pmatrix}\\begin{pmatrix} w_1\\\\ w_2\\\\ \\vdots\\\\ w_n\\end{pmatrix}\\\\\n",
+ "\\end{align}\n",
+ "\n",
+ "As can be seen from the expression above, not only is ${\\bf v}_{\\rm row}{\\bf w}_{\\rm \n",
+ "row}^{\\rm T}$ or ${\\bf v}_{\\rm col}^{\\rm T}{\\bf w}_{\\rm col}$ just as compact in notation as \n",
+ "${\\bf v}\\cdot{\\bf w}$, but it also offers the added benefit of explicitly specifying which of the two vectors\n",
+ "resides in \"row space\" (i.e., represented as a row vector) versus \"column space\" (i.e., represented as a column\n",
+ "vector). While this detail is seemingly inconsequential for our current definition of the dot product of two\n",
+ "real-valued vectors, drawing a distinction between row and column spaces is enormously important when working\n",
+ "with complex-valued vectors, or with the even more general entities with which quantum mechanics is built. While\n",
+ "the mathematical construction of quantum mechanics is beyond the scope of this lesson, we must nevertheless be\n",
+ "aware of the need to distinguish column and row space when building the software used to implement quantum\n",
+ "mechanics on a computer. From now onward, both for simplicity of notation and in order to maintain this \n",
+ "distinciton between row and column spaces, we will assume that all arbitrary vectors ${\\bf v}$ are column vectors,\n",
+ "and that their transposes ${\\bf v}^{\\rm T}$ are row vectors. The dot product is therefore assumed to be written\n",
+ "as\n",
+ "\n",
+ "$$\n",
+ "{\\bf v}\\cdot{\\bf w} = {\\bf v}^{\\rm T}{\\bf w}\n",
+ "$$\n",
+ "\n",
+ "#### Outer Product of Column/Row Vectors\n",
+ "\n",
+ "Considering the original expression we presented for the dot product,\n",
+ "\n",
+ "$$\n",
+ "{\\bf v}\\cdot{\\bf w} = \\sum_i v_i\\cdot w_i,\n",
+ "$$\n",
+ "\n",
+ "it should be clear that the dot product operation is commutative, i.e., ${\\bf v}\\cdot{\\bf w} =\n",
+ "{\\bf w}\\cdot{\\bf v}$. Now that we have rewritten the dot product as a matrix multiplication between the row vector\n",
+ "${\\bf v}^{\\rm T}$ and the column vector ${\\bf w}$, however, what would happen if we simply switched the order\n",
+ "of the two vectors in the product? In other words, if the dot product is given by ${\\bf v}^{\\rm T}{\\bf w}$, what\n",
+ "does the expression ${\\bf w}{\\bf v}^{\\rm T}$ yield?\n",
+ "\n",
+ "If the matrix product of a $1\\times N$ row vector and a $N\\times 1$ column vector yields a $1\\times 1$ matrix\n",
+ "(i.e., a scalar), then the matrix product of a $N\\times 1$ column vector and a $1\\times N$ row vector must\n",
+ "yield a $N\\times N$ matrix! This operation is called the _outer product,_ denoted with the $\\otimes$ symbol, and\n",
+ "is given by:\n",
+ "\n",
+ "$$\n",
+ "{\\bf v}\\otimes{\\bf w} = {\\bf v}{\\bf w}^{\\rm T} = \\begin{pmatrix} v_1 \\\\ v_2 \\\\ \\vdots \\\\ v_m\\end{pmatrix} \\begin{pmatrix}w_1 & w_2 & \\cdots & w_n\\end{pmatrix} = \\begin{pmatrix} v_1w_1 & v_1w_2 & \\cdots & v_1w_n\\\\\n",
+ "v_2w_1 & v_2w_2 & \\cdots & v_2w_n\\\\ \\vdots & \\vdots & \\ddots & \\vdots\\\\ v_mw_1 & v_mw_2 & \\cdots & v_mw_n\\end{pmatrix}\n",
+ "$$\n",
+ "\n",
+ "> Note: Just like we used a different notation within code to denote the inner product, we will denote the outer\n",
+ "product of two vectors, ${\\bf v}\\otimes{\\bf w}$, as `|v>\n",
+ "\n",
+ "\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Create np.array's for all matrices & vectors from above <==\n",
+ "# Define variables for numpy arrays of vectors v, w and matrices A, B from above\n",
+ "np_A = np.array(A)\n",
+ "np_B = np.array(B)\n",
+ "np_v = np.array(v)\n",
+ "np_w = np.array(w)\n",
+ "\n",
+ "# What types are these?\n",
+ "print(type(A))\n",
+ "print(type(np_A))\n",
+ "print(type(np_v))"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Array Information\n",
+ "\n",
+ "Unlike when we used basic Python `list`s to represent vectors and matrices above, NumPy arrays carry relevant\n",
+ "information, like their shape, around with them. So, instead of asking for `len(A)` and `len(A[0])` to determine\n",
+ "the shape of a matrix `A`, the shape of the NumPy array `np_A` is contained within the `np_A.shape` attribute.\n",
+ "Try it out below!"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 19,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "(3, 3)\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Array attributes are useful! <==\n",
+ "\n",
+ "print(np_A.shape)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Array Operations\n",
+ "\n",
+ "Unlike when using `list`-based representations of matrices and vectors, performing array operations is much more\n",
+ "straightforward with NumPy arrays because it is no longer necessary to operate on individual array elements. So,\n",
+ "instead of needing to iterate over each array element to perform, e.g., scalar multiplication, it is instead\n",
+ "possible to simply use the Python multiplication operator `*`, thanks to a useful NumPy trick called\n",
+ "_broadcasting_. \n",
+ "\n",
+ "In the cell below, evaluate the indicated expressions _without_ using `for` loops:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 20,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "[50 70 90]\n",
+ "[ 85 105 125]\n",
+ "[ 3.33333333 6.66666667 10. ]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Scalar array operations with NumPy Broadcasting <==\n",
+ "\n",
+ "# Evaluate v + w\n",
+ "tmp = np_v + np_w\n",
+ "print(tmp)\n",
+ "\n",
+ "# Evaluate 2 * w + 5\n",
+ "tmp = 2 * np_w + 5\n",
+ "print(tmp)\n",
+ "\n",
+ "# Evaluate v / 3\n",
+ "tmp = np_v / 3\n",
+ "print(tmp)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "In addition to the syntactic simplicity afforded by these NumPy-enabled array operations, they will also tend\n",
+ "to be much faster to execute than our by-hand solutions. To see this, we can use the Jupyter _magic function_\n",
+ "`%timeit`, which will report the time necessary to execute any line of code in a notebook. In order to make the\n",
+ "difference easier to see, as well, let's use a few large vectors which can be automatically generated by another\n",
+ "NumPy function, `np.random.random()`:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 21,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "By-hand r*v + s:\n",
+ "\t length-1000 vector:\n",
+ "730 µs ± 52.2 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)\n",
+ "\t length-10000 vector:\n",
+ "7.05 ms ± 70.7 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n",
+ "NumPy r*v + s:\n",
+ "\t length-1000 vector:\n",
+ "3.06 µs ± 79.6 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)\n",
+ "\t length-10000 vector:\n",
+ "10.3 µs ± 47.5 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Timing our vector operations vs NumPy <==\n",
+ "\n",
+ "# Define some big vectors\n",
+ "a = np.random.random(1000)\n",
+ "b = np.random.random(10000)\n",
+ "r = 3\n",
+ "s = 500\n",
+ "\n",
+ "print('By-hand r*v + s:')\n",
+ "print('\\t length-1000 vector:')\n",
+ "%timeit rvps(a, r, s)\n",
+ "print('\\t length-10000 vector:')\n",
+ "%timeit rvps(b, r, s)\n",
+ "\n",
+ "print('NumPy r*v + s:')\n",
+ "print('\\t length-1000 vector:')\n",
+ "%timeit r*a + s\n",
+ "print('\\t length-10000 vector:')\n",
+ "%timeit r*b + s\n"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "While the speeds of each of the above operations will change depending on your computer, the NumPy operations\n",
+ "should be ***much*** faster. As of writing this lesson, the NumPy operation of $r\\cdot{\\bf v} + s$ on my laptop\n",
+ "is approximately $200\\times$ faster than my by-hand solution for the length-1,000 vector, and approximately \n",
+ "$6,000\\times$ faster with the length-10,000 vector! \n",
+ "\n",
+ "As we will see, this difference in speed between a by-hand implementation and NumPy will become\n",
+ "even more drastic for the types of operations and objects which are used in molecular physics. But what exactly\n",
+ "are these objects?"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Tensors\n",
+ "\n",
+ "So far, we have seen that scalars, vectors, and matrices all share basic properties, and can interact with one\n",
+ "another through operations like addition and multiplication. Furthermore, we have seen that vectors and matrices\n",
+ "behave similarly, with their differences arising from the fact that vectors are one-dimensional arrays and\n",
+ "matrices are two-dimensional arrays. In a similar fashion, scalars could be considered to be 0-dimensional\n",
+ "arrays. So, if thus far we have considered the properties of 0-, 1-, and 2-dimensional arrays, what's to stop us\n",
+ "from extending our understanding to $N$-dimensional arrays? Before we do consider arbitrary-dimension arrays and\n",
+ "their properties, however, it is important to understand how and why they connect to the scalars, vectors, and\n",
+ "matrices we have already developed. Formally, the reason that scalars, vectors, and matrices all behave similarly\n",
+ "is because they are all examples of a more general type of object, which we will refer to as a _tensor_.\n",
+ "\n",
+ "_Tensors_ are a general class of mathematical entity related to vector spaces, which includes vectors and other\n",
+ "$N$-dimensional arrays, functions, and even operations like the derivative and the dot product. With this breadth\n",
+ "of different types of objects which are all technically tensors, we must have some way to denote tensors which\n",
+ "is broadly applicable. To this end, we will denote a tensor by using subscripted or superscripted indices:\n",
+ "\n",
+ "- Vectors, matrices, & $N$-dimensional arrays: $v_i$, $M_{ij}$, $T_{ij\\cdots k}$\n",
+ "- Functions, maps, & operators: $\\hat{f}_{xy}$, ${\\cal F}_i^j$, $\\hat{\\scr O}_{ij}$\n",
+ "\n",
+ "From our discussion of scalars, vectors, and matrices as 0-, 1-, & 2-dimensional arrays, we can see that the\n",
+ "\"dimension\" of the array is the same as the number of indices used to represent the tensor. To disambiguate\n",
+ "between the concept of array \"dimension\" and the dimension of a vector space (i.e., the number of basis vectors),\n",
+ "we will refer to the number of indices used to denote a tensor as its _rank_. So, scalars, vectors, and matrices\n",
+ "are rank-0, rank-1, and rank-2 tensors, respectively. With this new notation at our disposal, let's explore tensor\n",
+ "operations from the perspective of viewing tensors as multidimensional arrays. "
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Tensor Operations\n",
+ "\n",
+ "#### Elementwise Tensor Operations\n",
+ "\n",
+ "Just like for vectors and matrices, rank-$N$ tensors also have defined a scalar multiplication and addition\n",
+ "operations, as well as elementwise array operations. To explore this, use the cell below to first define two\n",
+ "rank-4 NumPy arrays $M_{pqrs}$ and $N_{pqrs}$, before evaluating the indicated expressions.\n",
+ "\n",
+ "> Note: When two tensors use the same indices, they are assumed to have identical shape.\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 22,
+ "metadata": {
+ "scrolled": true
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "5*M + 2 =\n",
+ "[[[[4.41647375 3.75037663]\n",
+ " [3.14386216 5.75129335]]\n",
+ "\n",
+ " [[2.77291831 5.17161798]\n",
+ " [3.86402171 5.45887502]]]\n",
+ "\n",
+ "\n",
+ " [[[2.5422852 6.3967675 ]\n",
+ " [5.47417575 3.66519298]]\n",
+ "\n",
+ " [[2.92498233 5.29089265]\n",
+ " [3.01843885 5.67434326]]]]\n",
+ "\n",
+ "M + N =\n",
+ "[[[[0.83138173 0.87362918]\n",
+ " [1.07625849 1.60398442]]\n",
+ "\n",
+ " [[0.91804881 1.54211613]\n",
+ " [1.21460153 1.65162547]]]\n",
+ "\n",
+ "\n",
+ " [[[0.45053935 0.89339571]\n",
+ " [0.91410228 0.90504179]]\n",
+ "\n",
+ " [[1.13398664 1.56466725]\n",
+ " [1.0049532 1.08667003]]]]\n",
+ "\n",
+ "M*N + 10 =\n",
+ "[[[[10.16822861 10.18328328]\n",
+ " [10.19388145 10.64051515]]\n",
+ "\n",
+ " [[10.11801924 10.57583422]\n",
+ " [10.31382565 10.66400056]]]\n",
+ "\n",
+ "\n",
+ " [[[10.03710123 10.01234806]\n",
+ " [10.15235451 10.19049914]]\n",
+ "\n",
+ " [[10.17555983 10.59663141]\n",
+ " [10.16320797 10.2585278 ]]]]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Scalar & Elementwise Tensor Operations <==\n",
+ "# Declare two rank-4 tensors, Mpqrs & Npqrs, using np.random.random()\n",
+ "Mpqrs = np.random.random((2,2,2,2)) # Limit the dimensions to be no more than length-3 each\n",
+ "Npqrs = np.random.random((2,2,2,2)) # Make sure the dimensions are the same as Mpqrs!\n",
+ "\n",
+ "# Evaluate 5 * M + 2\n",
+ "tmp = 5 * Mpqrs + 2\n",
+ "print(\"5*M + 2 =\")\n",
+ "print(tmp)\n",
+ "\n",
+ "# Evaluate M + N\n",
+ "tmp = Mpqrs + Npqrs\n",
+ "print(\"\\nM + N =\")\n",
+ "print(tmp)\n",
+ "\n",
+ "# Evaluate M*N + 10; recall `*` indicates the elementwise product\n",
+ "tmp = Mpqrs * Npqrs + 10\n",
+ "print(\"\\nM*N + 10 =\")\n",
+ "print(tmp)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "#### Tensor Contractions: Generalized Array Multiplication\n",
+ "\n",
+ "The way we defined the matrix product above, we may only multiply compatible arrays which share the same \"inner\"\n",
+ "dimensions to yield a matrix with the \"outer\" dimensions. Let's take a closer look at the triple-summation\n",
+ "form of the matrix multiplication operation:\n",
+ "\n",
+ "$${\\bf C} = {\\bf A}\\times{\\bf B} = \\sum_{i=1}^{M}\\sum_{k=1}^{N}\\sum_{j=1}^{P}A_{ik}B_{kj} = C_{ij}$$\n",
+ "\n",
+ "Here, we can see that we are multiplying two rank-2 tensors, $A_{ik}$ and $B_{kj}$, to produce another rank-2\n",
+ "tensor, $C_{ij}$. This and other non-elementwise tensor-tensor multiplications will be referred to as _tensor\n",
+ "contractions,_ which can be thought of as generalized matrix-matrix multiplications which occur over particular\n",
+ "tensor indices. Thanks to the fact that we denote tensors based on their indices, however, we no longer need\n",
+ "to explicitly concern ourselves with the _order_ of the indices in a contraction; for example, the contraction\n",
+ "above could therefore be rewritten as\n",
+ "\n",
+ "$$C_{ij} = \\sum_{i=1}^{M}\\sum_{k=1}^{N}\\sum_{j=1}^{P}A_{ik}B_{jk} = {\\bf A}\\times {\\bf B}^{\\rm T},$$\n",
+ "\n",
+ "which may not be allowed by the shapes of the matrices ${\\bf A}$ and ${\\bf B}$ (if, e.g., ${\\bf A}$ is\n",
+ "$3\\times 4$ but ${\\bf B}$ is $4\\times 3$). Clearly, writing even a basic matrix multiplication\n",
+ "as a sum over common indices offers increased flexibility over the conventional definition of matrix\n",
+ "multiplication. By examining these summation expressions further, it should be apparent that only terms where\n",
+ "values of the index $k$ are shared contribute to the summation. Therefore, it is acceptable to remove the\n",
+ "explicit summations over indices $i$ and $j$, instead only retaining the summation over $k$:\n",
+ "\n",
+ "$$C_{ij} = \\sum_{i=1}^{M}\\sum_{k=1}^{N}\\sum_{j=1}^{P}A_{ik}B_{jk} = \\sum_{k} A_{ik}B_{kj}.$$\n",
+ "\n",
+ "Because it is understood, however, that only the terms involving shared values for the index $k$ are retained in\n",
+ "the summation, it is also convenient _not_ to write the sum at all:\n",
+ "\n",
+ "$$C_{ij} = \\sum_{i=1}^{M}\\sum_{k=1}^{N}\\sum_{j=1}^{P}A_{ik}B_{jk} = \\sum_{k} A_{ik}B_{kj} = A_{ik}B_{kj}.$$\n",
+ "\n",
+ "In this step, we have leveraged the _Einstein summation convention_ to simplify the notation for our tensor\n",
+ "contraction:\n",
+ "\n",
+ "> Einstein summation convention: In a tensor expression, repeated indices are assumed to be summed over.\n",
+ "\n",
+ "Let's use this convention to redefine the array multiplications we introduced above in tensor notation!\n",
+ "\n",
+ "| Product Type | Array Notation | Einstein Summation | Example Shape |\n",
+ "|--------------|----------------|--------------------|---------------|\n",
+ "| Vector Inner Product | ${\\bf v}\\cdot{\\bf w}$ | $v_i w_i\\rightarrow r$ | (1, $N$) x (1, $N$) $\\rightarrow$ (1, 1) |\n",
+ "| Vector Outer Product | ${\\bf v}\\otimes{\\bf w}$ | $v_i w_j\\rightarrow M_{ij}$ | (1, $N$) x (1, $M$) $\\rightarrow$ ($N$, $M$) |\n",
+ "| Matrix Inner Product | ${\\bf A}\\cdot{\\bf B}^{\\rm T}$ | $A_{ik}B_{kj}\\rightarrow C_{ij}$ | (2, 3) x (3, 4) $\\rightarrow$ (2, 4) |\n",
+ "| Matrix Outer Product | ${\\bf A}\\cdot{\\bf B}$ | $A_{ij}B_{ik}\\rightarrow C_{jk}$ | (2, 8) x (2, 5) $\\rightarrow$ (8, 5) |\n",
+ "\n",
+ "> Note: For the vector outer product, no index is shared between the two rank-1 tensors ${\\bf v}$ and ${\\bf w}$;\n",
+ "therefore, the resulting array is of shape $N\\times M$. \n",
+ "\n",
+ "While each of these product types are simple to perform using standard matrix-vector or matrix-matrix \n",
+ "multiplication (i.e., by using `np.dot()`), it is challenging to do so for more complex tensor contractions. \n",
+ "Instead, NumPy contains a function which allows for arbitrary contractions according to Einstein summation\n",
+ "convention, `np.einsum()`, which takes a \"map\" of the indices involved in the contraction as an argument:\n",
+ "\n",
+ "| Product Type | Array Notation | Einstein Summation | Example Shape | `np.einsum()` Call |\n",
+ "|--------------|----------------|--------------------|---------------|--------------------|\n",
+ "| Vector Inner Product | ${\\bf v}\\cdot{\\bf w}$ | $v_i w_i\\rightarrow r$ | (1, $N$) x (1, $N$) $\\rightarrow$ (1, 1) | `np.einsum('i,i->', v, w)` |\n",
+ "| Vector Outer Product | ${\\bf v}\\otimes{\\bf w}$ | $v_i w_j\\rightarrow M_{ij}$ | (1, $N$) x (1, $M$) $\\rightarrow$ ($N$, $M$) | `np.einsum('i,j->ij', v, w)` |\n",
+ "| Matrix Inner Product | ${\\bf A}\\cdot{\\bf B}^{\\rm T}$ | $A_{ik}B_{kj}\\rightarrow C_{ij}$ | (2, 3) x (3, 4) $\\rightarrow$ (2, 4) | `np.einsum('ik,kj->ij', A, B)` |\n",
+ "| Matrix Outer Product | ${\\bf A}\\cdot{\\bf B}$ | $A_{ij}B_{ik}\\rightarrow C_{jk}$ | (2, 8) x (2, 5) $\\rightarrow$ (8, 5) | `np.einsum('ij,ik->jk', A, B)` |\n",
+ "\n",
+ "In the cell below, use `np.einsum` to evaluate the indicated tensor expression:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 23,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Aij Babjd -> [[[[1.98573675 0.9768675 ]\n",
+ " [1.67181843 1.24962792]\n",
+ " [1.31219973 1.40674247]\n",
+ " [1.53483055 1.13157685]]\n",
+ "\n",
+ " [[2.07156415 0.71159499]\n",
+ " [1.74408436 1.77837014]\n",
+ " [1.95023464 1.30448452]\n",
+ " [1.57524559 1.25462834]]]\n",
+ "\n",
+ "\n",
+ " [[[1.89545155 0.56626646]\n",
+ " [1.7433513 1.23860222]\n",
+ " [1.34911293 1.19384743]\n",
+ " [1.20446862 1.2198087 ]]\n",
+ "\n",
+ " [[1.96578292 0.85914056]\n",
+ " [1.72098592 1.93628759]\n",
+ " [1.78421445 1.6043322 ]\n",
+ " [1.28763136 1.41790362]]]\n",
+ "\n",
+ "\n",
+ " [[[1.01268171 0.61542484]\n",
+ " [0.86207709 0.60981275]\n",
+ " [0.60731096 0.90640351]\n",
+ " [0.74773556 0.60384708]]\n",
+ "\n",
+ " [[1.17139199 0.38207571]\n",
+ " [0.7094918 0.96227062]\n",
+ " [0.93628894 0.75110666]\n",
+ " [0.90975923 0.60766486]]]]\n",
+ "Aij Cii -> [0.80053656 1.06485495 0.59884852 0.7351457 0.57408373]\n",
+ " = 3.4502643890304223\n",
+ "Cij Aik -> [[0.62823687 0.45523464 0.25916895 0.42418631 0.38007634]\n",
+ " [0.51767451 0.92268938 0.50463239 0.57300983 0.41893123]\n",
+ " [0.72852992 0.55776695 0.28036196 0.51473638 0.45923349]]\n",
+ "|w> [[0.69907343 1.0574029 1.03984393 1.53180095 0.8408721 ]\n",
+ " [0.90707283 1.67779216 1.51702816 1.06498009 1.66150868]\n",
+ " [0.43546324 1.20840504 1.33397953 1.31771667 1.08074436]\n",
+ " [0.80503167 1.31178364 0.6207049 0.78337707 1.44639185]]\n",
+ "Bijkl -> Bikjl: [[[[0.65580342 0.80864031]\n",
+ " [0.33742978 0.40476387]\n",
+ " [0.25037273 0.8655413 ]\n",
+ " [0.69292722 0.27258856]]\n",
+ "\n",
+ " [[0.85221854 0.00511214]\n",
+ " [0.70358522 0.55241835]\n",
+ " [0.64338638 0.27917535]\n",
+ " [0.57668644 0.52283762]]\n",
+ "\n",
+ " [[0.55139732 0.12695251]\n",
+ " [0.72221488 0.28717171]\n",
+ " [0.34391793 0.64108926]\n",
+ " [0.05445456 0.54024296]]\n",
+ "\n",
+ " [[0.7835582 0.54669115]\n",
+ " [0.65827243 0.91237672]\n",
+ " [0.73718189 0.41130809]\n",
+ " [0.60731869 0.38483208]]\n",
+ "\n",
+ " [[0.75348413 0.0344589 ]\n",
+ " [0.90368564 0.01317345]\n",
+ " [0.4085451 0.50336814]\n",
+ " [0.46735079 0.56421452]]]\n",
+ "\n",
+ "\n",
+ " [[[0.77209139 0.04327001]\n",
+ " [0.1730143 0.56964305]\n",
+ " [0.63463836 0.1850905 ]\n",
+ " [0.73131845 0.11210445]]\n",
+ "\n",
+ " [[0.60975713 0.20518435]\n",
+ " [0.96172384 0.97420694]\n",
+ " [0.85398143 0.56501866]\n",
+ " [0.2508107 0.7350972 ]]\n",
+ "\n",
+ " [[0.8810984 0.48844661]\n",
+ " [0.22162479 0.79481328]\n",
+ " [0.32368201 0.99006159]\n",
+ " [0.53085927 0.56625034]]\n",
+ "\n",
+ " [[0.94543437 0.74824275]\n",
+ " [0.71779684 0.40670766]\n",
+ " [0.89988935 0.58053478]\n",
+ " [0.80381809 0.17605838]]\n",
+ "\n",
+ " [[0.80844891 0.08738797]\n",
+ " [0.94404491 0.75782303]\n",
+ " [0.61282246 0.67219926]\n",
+ " [0.58341149 0.97904106]]]]\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Practice with Einsum <==\n",
+ "# Declaring some tensors of various shape\n",
+ "A = np.random.random((3,5))\n",
+ "B = np.random.random((2,4,5,2))\n",
+ "C = np.random.random((3,3))\n",
+ "v = np.random.random((10,))\n",
+ "w = np.random.random((10,))\n",
+ "\n",
+ "# Tensor contraction A_{ij}B_{abjd}\n",
+ "AijBabjd = np.einsum('ij,abjd->iabd', A, B)\n",
+ "print(f\"Aij Babjd -> {AijBabjd}\")\n",
+ "\n",
+ "# Tensor contraction A_{ij}C_{ii}\n",
+ "AijCii = np.einsum('ij,ii->j', A, C)\n",
+ "print(f\"Aij Cii -> {AijCii}\")\n",
+ "\n",
+ "# Inner product \n",
+ "vdotw = np.einsum('i,i->', v, w)\n",
+ "print(f\" = {vdotw}\")\n",
+ "\n",
+ "# Tensor contraction C_{ij}A_{ik}\n",
+ "CijAik = np.einsum('ij,ik->jk', C, A)\n",
+ "print(f\"Cij Aik -> {CijAik}\")\n",
+ "\n",
+ "# Outer product of w with v\n",
+ "wouterv = np.einsum('i,j->ij', w, v)\n",
+ "print(f\"|w> B_{jk}\n",
+ "Bjk = np.einsum('ijki->jk', B)\n",
+ "print(f\"Bijki -> {Bjk}\")\n",
+ "\n",
+ "# Transpose operation B_{ijkl}->B_{ikjl}\n",
+ "Bikjl = np.einsum('ijkl->ikjl', B)\n",
+ "print(f\"Bijkl -> Bikjl: {Bikjl}\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "The single biggest benefit to `np.einsum` is how explicitly the contractions are represented. As an example,\n",
+ "consider the following contractions between a rank-4 tensor $I$ and a rank-2 tensor $D$:\n",
+ "$$J_{pq} = I_{pqrs}D_{rs}$$\n",
+ "$$K_{pq} = I_{prqs}D_{rs}$$\n",
+ "\n",
+ "While it is not obvious how to perform these contractions with `np.dot()`, these operations are simple to\n",
+ "translate into calls to `np.einsum()`. In the cell below, try it out:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 24,
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "(12, 12) (12, 12)\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Another example of Einsum simplicity <==\n",
+ "I = np.random.random((12, 12, 12, 12))\n",
+ "D = np.random.random((12,12))\n",
+ "\n",
+ "# Use einsum to compute J and K using the expressions above. \n",
+ "# Make sure you pay attention to the index ordering!\n",
+ "\n",
+ "J = np.einsum('pqrs,rs->pq', I, D)\n",
+ "K = np.einsum('prqs,rs->pq', I, D)\n",
+ "print(J.shape, K.shape)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Computational Efficiency of Tensor Contraction Engines"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 25,
+ "metadata": {
+ "scrolled": true
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Timing our matrix multiply:\n",
+ "CPU times: user 1min 49s, sys: 350 ms, total: 1min 49s\n",
+ "Wall time: 1min 50s\n",
+ "Timing np.einsum:\n",
+ "CPU times: user 42.7 ms, sys: 844 µs, total: 43.5 ms\n",
+ "Wall time: 43.2 ms\n",
+ "Timing np.dot:\n",
+ "CPU times: user 14.4 ms, sys: 2.49 ms, total: 16.9 ms\n",
+ "Wall time: 9.69 ms\n"
+ ]
+ }
+ ],
+ "source": [
+ "# ==> Declare large-ish matrices for timings <==\n",
+ "\n",
+ "A = np.random.random((500,500))\n",
+ "B = np.random.random((500,500))\n",
+ "\n",
+ "# Our hand-written matrix multiply\n",
+ "print('Timing our matrix multiply:')\n",
+ "%time mm_C = MM(A, B)\n",
+ "\n",
+ "# Einsum\n",
+ "print('Timing np.einsum:')\n",
+ "%time es_C = np.einsum('ik,kj->ij', A, B)\n",
+ "\n",
+ "# Dot product\n",
+ "print('Timing np.dot:')\n",
+ "%time dot_C = A.dot(B)\n"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Comparing Contraction Engines"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "##### Student Answer Box\n",
+ "1. Based on the experiences and use cases above, order the three contraction engines based on the following\n",
+ "factors from \"best\" to \"worst\". Justify your orderings.\n",
+ " 1. Computational efficiency (speed)\n",
+ " - `np.dot` > `np.einsum` >>> manual Python loops\n",
+ " - The timings are pretty clear.\n",
+ " 2. Code clarity & readability\n",
+ " - `np.einsum` > `np.dot` $\\sim$ manual Python loops\n",
+ " 3. Engine flexibility\n",
+ " - `np.einsum` > `np.dot` >>> manual Python loops\n",
+ " \n",
+ "2. Based on your orderings, recommend a use case for each contraction engine. Justify your recommendation.\n",
+ " 1. Manual Python loops\n",
+ " - Either don't use or only use to teach the matrix multiplication formula.\n",
+ " 2. NumPy Einsum\n",
+ " - Complicated contraction, etc. and want to be explicit while maintaining decent efficiency\n",
+ " 3. NumPy Dot\n",
+ " - Any time that readability is not as important as speed."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": []
+ }
+ ],
+ "metadata": {
+ "kernelspec": {
+ "display_name": "Python 3",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.8.5"
+ },
+ "latex_envs": {
+ "LaTeX_envs_menu_present": true,
+ "autoclose": false,
+ "autocomplete": true,
+ "bibliofile": "biblio.bib",
+ "cite_by": "apalike",
+ "current_citInitial": 1,
+ "eqLabelWithNumbers": true,
+ "eqNumInitial": 1,
+ "hotkeys": {
+ "equation": "Ctrl-E",
+ "itemize": "Ctrl-I"
+ },
+ "labels_anchors": false,
+ "latex_user_defs": false,
+ "report_style_numbering": true,
+ "user_envs_cfg": true
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}
diff --git a/Tutorials/02_Linear_Algebra/README.md b/Tutorials/02_Linear_Algebra/README.md
index 9f1c7238..bb2e8715 100644
--- a/Tutorials/02_Linear_Algebra/README.md
+++ b/Tutorials/02_Linear_Algebra/README.md
@@ -1,4 +1,16 @@
-## Linear Algebra
+A Brief Introduction to Linear Algebra
+======================================
-**placeholder**
+This module seeks to provide a brief overview of the linear algebra topics
+which are relevant to the implementation of quantum chemistry methods in Python
+that we will begin in Module 03, with the following tutorials:
+
+- (2a) Vectors, Matrices, & Tensors and their Operations: Introduces vectors,
+matrices, and tensors, as well as their representations in Python, and
+discusses the basics of their operations.
+
+
+#### Planned Tutorials:
+
+- Eigenproblems & Operator Theory