Dr. Robert van de Geijn: The title of this course is Linear Algebra
Foundations To Frontiers.
What we would like to do is already give you a glimpse of the frontier in
this first week.
The frontier we're going to show you is how we, as part of our Flame
project here at the University of Texas, encode algorithms in code.
Now our philosophy is that the algorithm should be expressed in code
as much as possible in a way that reflects how one
thinks about the problem.
So in the last unit we showed how to think about the algorithm in terms of
sub vectors What we're now going to show you is how to represent that in
code in terms of sub vectors.
Let's have a look.
Here is the dashboard for your Python notebooks, and what I would like for
you to do is open up 1.6.3 Programming without indices for dot product.
Now the first thing we always do is to make a copy.
So let me go ahead and make a copy.
All right, now in this notebook if you read through it, the first thing we do
is we create factors x and y as matrices.
And then we reiterate the simple routine that we wrote in an earlier
unit, that computes the dot product of x and y outputting
the result as alpha.
So let's execute that box, and then let's just test that this actually
does what it's supposed to do.
What we're going to introduce you to next is how we now encode this in a
way that very closely resembles this algorithm as we expressed
it in the last unit.
Let's put this algorithm on the side here, so that we can
look at it as we code.
And let me just move this over a little bit, and now we have the
algorithm on the left.
And what we would like to do is create code that very closely resembles this
algorithm on the left, but in Python code.
Now how are we going to do this.
OK this is where the magic starts.
What we're going to do is visit the web page that we call
the spark web page.
And that web page is going to help us code the algorithm on the left in a
way that closely mirrors an algorithm, but in Python.
So let's click on this.
And what we have now is a menu.
And we're going to focus first on the part on the left.
And I strongly suggest that you follow along as we do this.
Let me make this a little bigger so that you can see what's going on.
The first thing we need to do is give the routine that we're
going to write a name.
So we're going to call it dot, and let's make sure that we use lowercase
letters here.
The type that we're going to implement is called unblocked.
Now what unblocked versus blocked means we'll see later.
So let's just leave this as unblocked.
The variant number has to do with the fact that there are often many
different algorithms that compute the same operation.
So we could've numbered these, for the moment we're going to
leave this as none.
Now how many different operands are involved in the dot product?
Well there's the vector x, there's the vector y, and there's the scalar alpha
that is the result.
So that means there are three.
And what we get is we get a menu for the different operands that are going
to be passed in.
Now the first operand we would like to be x.
The second operand, y.
And the third operand is going to be alpha, which in this
menu we pick as alpha.
Now x is a vector, y is a vector, and a is a scalar.
And that's going to make it show up as alpha.
X and y are merely inputs, we're not going to change those vectors.
And alpha, or a, is going to be both input and output because
that's how we work.
The output line, which we pick as FlamePy, which is an application
programming interface that we have created on top of Python.
What this exactly means will become clear as we move on.
Now if you want you can go and put your name in here and your email.
I'll just put in my name because I really don't you people to email me.
There are just too many of you, I'm sorry.
And then comes a very important part, which has to do with how we've marched
through the vector.
Now if we focus on the algorithm on the left here, what we notice is that
x is being partitioned into a top part and a bottom part, and the same for y.
When we start, x top has no components, 0 elements, and when we're
finished x top consists of the entire vector.
So what this means to us is that we march through this vector from top to
bottom as we compute.
So over here, for vectors x and y we choose to go from top to bottom, t
arrow b means top to bottom, and the same for y.
So as we were filling out this menu on the right, we ended up with an
implementation in Python.
Now if you take that implementation and you put it next to the algorithm,
you notice some similarities.
Let's see if we can move this all the way here, and we could put these next
to each other.
In the algorithm we start by saying that x needs to be partitioned into a
top part and a bottom part, and the same for y, where x top and y top have
no components initially.
Over here we do the same thing with a subroutine call where we have the
flame package.
And the flame package has a routine that takes a vector x and partitions
it into a top part and a bottom part, where this here indicates that the top
part of x initially has zero components.
That means that x, the bottom part of x is all of x.
Same thing for x.
In the algorithm we next check whether x top has as many components as x,
because then we're finished.
On the right we do the same thing by extracting the shape, the row and
column dimension of x top and x.
And comparing the row dimension, and to see whether we're done.
Next we have a routine that takes the top and bottom part and peel one
element off of x bottom, and re-names the sub vectors x0, chil1, and x2.
And the same thing for y.
That corresponds to this re-partition that we have on the left.
At the bottom of the algorithm, we say continue by redefining x top to
consist of x0 and chil1.
And the same thing for y top, y0 and psil1.
On the right we have a retain that says take the middle element, chi 1,
add it to the top and call that x top.
x2, becomes x bottom.
So we have a way of implementing this continue with that
we have in our algorithm.
OK so what do we do?
We highlight all of this code, we copy it, we take it over here, we insert it
into this box where it says insert code here, and that are a couple of
things that we need to do.
Over here it also says that alpha must be set to 0.
So what we need to do over here is set alpha equal to 0.
Except that we are passing alpha in as a matrix, so what we really should do
is set the 0,0 element of that matrix to 0.
And then that web page doesn't really know what operation
you're trying to implement.
So right here what we need to do is look and see what update
must be done to alpha.
And then we need to insert that right here.
So we're going to say, alpha of 0,0 because we want to actually override
that element in the matrix, is equal to chi1 times psi1, plus alpha 0,0.
Now we need to execute this block, this cell.
So we do that.
And then here we have a cell in which we have code that tests to see whether
we get the right answer.
What's happening here, alpha is first declared as a matrix that's 1 by 1,
with the element minus 1 in it, just to give it some value.
It then calls the dot routine which should override alpha with the dot
product of x and y.
And then it compares here the result to what we would compute if we did x
transpose times y.
So let's execute that cell.
And what we find is that we get the right answer.
Welcome to the wonderful world of how we program in the FLAME group.
What we actually do is we code in C, and we've actually written libraries
that have hundreds of thousands of lines of code this way.
These codes are highly reliable because it's very difficult to make
errors the way we code, and they actually tend to also achieve very
high performance.
And how they achieve high performance we'll talk about later in this course.