The following content is provided under a Creative
Commons license.
Your support will help MIT OpenCourseWare continue to
offer high quality educational resources for free.
To make a donation or view additional materials from
hundreds of MIT courses, visit MIT OpenCourseWare at
ocw.mit.edu.
PROFESSOR: Good, we're going to take a detour today into
the realm of algorithms.
So when you're trying to make code go fast, of course,
there's no holds barred.
You can use whatever you need to in order
to make it go fast.
Today we're going to talk a little bit in a more
principled way about the memory hierarchy.
And to do that we're going to introduce what we call the
ideal-cache model.
So as you know most caches are hacked together to try to
provide something that will cache well while still making
it easy to build and fast to build.
The ideal-cache model is a pretty nice
beast if we had them.
It's got a two-level hierarchy.
It's got a cache that has m bytes that are organized in to
b byte cache-lines.
So each block is b bytes.
And it's fully associative.
So you recall, that means that any line can
go anywhere in cache.
And probably the most impressive aspect of an
ideal-cache is that it has an optimal
omniscient replacement algorithm.
So what it does, is it figures out when it needs to kick
something out of cache, it says, what is the absolutely
best thing you could possibly kick out of cache.
And it does that one.
Looking into the future if need be.
It says, oh is this going to be accessed
a lot in the future?
I think I'll keep this one.
Let's throw out this one.
I know it's never going to be used again.
So it has that omniscient character to it.
The performance measures we're going to look at in this
model, the first one is what we call the work.
And that's just the ordinary serial running time if you
just ran the code on one processor and counted up
essentially how many processor instructions you would do.
That's essentially the work.
The second measure, which is the one that's much more
interesting, is cache misses.
So the work has to do with the processor.
The cache misses has to do with what moves between these
two levels of memory.
So in this case, what we're interested in doing is, how
often do I try to access something.
It's not in the cache.
I have to go to main memory and bring it back into cache.
And so that's what we'll be counting in this model.
So it's reasonable to ask how reasonable ideal caches are.
In particular the assumption of omniscient replacement,
that's pretty powerful stuff.
Well it turns out there's a great lemma due to Slater and
Tarjan that says essentially the following.
Suppose that you have an algorithm that incurs q cache
misses on an ideal cache of size n.
So you ran the algorithm on your machine, you had a
cache of size n.
Then, if instead of having an ideal cache, you have a fully
associative cache of size two m and use the least recently
used replacement policy.
So you always, whenever you're kicking something out of
cache, you kick out the thing that has been touched the
longest ago in the past.
Then it incurs at most 2Q cache misses.
So what that says is that LRU is to with it constant factors
essentially the same as optimal.
Really quite a remarkable result.
Who's taking 6046?
You've just seen this, right?
Yeah, OK.
Just seen this result in 6046.
See, I do talk to my colleagues occasionally.
So then something about how this is proved.
And what's important here is that really it just says, OK,
Yeah you could dither on the constants, but basically
whether you choose LRU or choose ideal cache with the
omniscient replacement, asymptotically you're not
going to be off at all.
So for most asymptotic analyses, you can assume
optimal or LRU replacement as convenient.
And the typical way that you do convenience is if you're
looking at upper bounds.
So you're trying to show that a particular algorithm is
good, then what you do is you assume optimal replacement.
If you're trying to show that some algorithm is bad, then
what you do is assume that it's LRU to get a lower bound.
Because then you can reason more easily about what's
actually in memory.
Because you just say, oh we'll just keep the least
recently used one.
So you tend to use the two for upper bounds and lower bounds.
Now, the way this relates to software
engineering is as follows.
If you're developing a really fast algorithm, it's going to
start from a theoretically sound algorithm.
And from that then you have to engineer for detailed
performance.
So you have to take into account things like real
caches are not fully associative.
That loads and stores, for example, have different cost
with respect to bandwidth and latency.
So whether you miss some a load or miss on a store,
there's a different impact.
But these are all the tuning.
And as you know, those constant factors can sometimes
add up to dramatic numbers, orders of magnitude.
And so it's important to do that software engineering.
But starting from a good theoretical basis means that
you actually have an algorithm that is going to work well
across a large variety of real situations.
Now, there's one other assumption we tend to make
when we're dealing with ideal caches, and that's called the
tall-cache assumption.
So what the tall-cache assumption says, is that I you
have at least as many lines of cache, essentially, in your
cache, as you have bytes in the line.
So it says the cache is tall.
In other words, this dimension here is bigger than this
dimension here.
And in particular, you want that to be true for where we
have some constant here of slop that we can throw in.
Yes, question.
AUDIENCE: Does that [INAUDIBLE]
associatively make the cache shorter here.
PROFESSOR: Yes, so this is basically assuming
everything is ideal.
We're going to go back.
When you engineer things, you have to deal with the fact
that things aren't ideal.
But usually that's just a little bit of a tweak on the
actual ideal algorithm.
And for many programs, the ideal algorithm you don't
actually have to tweak at all to get a
good practical algorithm.
So here is just saying the cache should be tall.
Now, just as an example, if we look at the machines that
we're using, the cache-line length is 64 bytes.
The L1 cache size is 32 kilobytes.
And of course, for L1 it's is 32 kilobytes and for L2 and
L3, it's even bigger.
It's even taller.
Because they also have 64k line.
So this is a fairly reasonable assumption to make, that you
have more lines in your cache then essentially the length,
the number of items you can put on a cache line.
Now why is this an important assumption?
So what's wrong with short caches?
So we're going to look at, surprise, surprise, matrix
multiplication.
Which, by the end of this class you will learn more
algorithms than matrix multiplication.
But it is a good one to illustrate things.
So the idea here is, suppose that you have an
n by n matrix here.
And you don't have this tall-cache assumption.
So where your cache is short.
You have a lot of bytes in a line, but very few lines.
Then even if the size of your matrix fits, in
principle, in the cache.
In other words, n squared is less than m by more than a
constant amount.
So you'd say, Oh gee, that ought to fit in.
If you have a short cache it doesn't necessarily fit in
because your length n here is going to be shorter than the
number of bytes on a line.
However, if you have a tall cache, it's always the case
that if the matrix size is smaller than the cache size by
a certain amount, then the matrix will fit in the cache.
OK, question?
AUDIENCE: Why wouldn't you fit more than one
row per cache line?
PROFESSOR: Well the issue is you may not have control over
the way this is laid out.
So, for example, if this is row-major order, and this is a
submatrix of a much bigger matrix, then you may not have
the freedom to be using these.
But if you have the tall-cache assumption, then any section
you pull out is going to fit.
As long as the data fits mathematically in the cache,
it will fit practically in the cache if you have the
tall-cache assumption.
Whereas if it's short, you basically end up with the
cache lines being long and you not having any flexibility as
to where the data goes.
So this is sort of a--
So any questions about that before we get
into the use of this?
We're going to see the use of this and where it comes up.
So one of the things is that, if it does fit in, then it
takes, at most, size of the matrix divided by the cache
line size misses to load it in.
So this is linear time in the cache world.
Linear time says, you should only take one cache fault for
every line of cache.
And so that's what you'll have here if you have a tall cache.
You'll have n squared over b cache misses to load in n
square data.
And that's good.
OK, good.
So let's take on the problem of multiplying matrices.
We're going to look at square matrices because they're
easier to think about than rectangular ones.
But almost everything I say today will relate to
rectangular matrices as well.
And it we'll generalize beyond matrices as
we'll see next time.
So here's a typical code for multiplying matrices.
It's not the most efficient code in the world, but it's
good enough to illustrate what I want to show you.
So the first thing is, what is the work of this algorithm?
This is, by the way, the softball question.
What's the work?
So the work, remember, is just if you're analyzing it just
like processor forget about caches and so forth.
AUDIENCE: n cubed.
PROFESSOR: n cubed, right.
Because there's a triply nested loop going up to n and
you're doing constant work in the middle.
So it's n times n times 1.
n cubed work.
That was easy.
Now let's analyze caches.
So we're going to look at row major.
I'm only going to illustrate the cache lines on this side
because B is where all the action is.
So we're going to analyze two cases when the matrix doesn't
fit in the cache.
If the matrix fits in the cache, then there's nothing to
analyze, at some level.
So we're going to look at the cases where the matrix doesn't
fit in the cache.
And the first one is going to be when the side of the matrix
is bigger than m over b.
So remember, m over b is the height of our cache, the
number of lines in our cache.
So let's assume for this, now I have a choice of assuming
optimal omniscient replacement or assuming LRU.
Since I want to show this is bad, I'm going to assume LRU.
Could somebody please close the back door there?
Because it's we're getting some noise in from out there.
Thank you.
So let's assume LRU.
So what happens in the code is basically I go across a row of
A, while I go down a column of B. And now if I'm using LRU,
what's happening?
I read in this cache block and this one,
then this one et cetera.
And if n is bigger than M/B and I'm using least recently
used, by the time I get down to the bottom here, what's
happened the first cache line?
First cache block?
It's out of there.
It's out of there if I used LRU.
So therefore, what happens is I took a miss on
every one of those.
And then when I go to the second one, I take a miss on
every one again.
And so as I keep going through, every access to B
causes a miss throughout the whole accessing of B. Now go
on to the second row A and I had the same thing repeats.
So therefore, the number of cache misses is order n cubed
since we miss on matrix B on every access.
OK, question.
AUDIENCE: I know that you said it's e.
Does B push out due to conflict
misses or capacity misses?
PROFESSOR: So in this case they're capacity misses that
we're talking about here.
So there's no conflict misses in a fully associative cache.
Conflict misses occurs because of direct mapping.
So there's no conflict misses in what I'm going to be
talking about today.
That is an extra concern that you have for real caches, not
a concern when you have a fully associative cache.
AUDIENCE: [INAUDIBLE]
n needs to be bigger than B?
PROFESSOR: So the number of lines to fit in my cache is m
over b, right?
AUDIENCE: So can't you put multiple units of data--
PROFESSOR: Well there are multiple units of data.
But notice this is row major, what's stored here
is B11, B12, B13.
That's stored here.
The way I'm going through the access, I'm going down the
columns of B. So by the time to get up to the top again,
that cache block is no longer there.
And so when I access B12, assuming indexing from one or
whatever, this block is no longer in cache.
Because LRU would say, somewhere along here I hit the
limit of my size of cache, let's say around here.
Then when this one goes in, that one goes out.
When the next one goes in, the next one goes out et cetera
using the least recently used replacement.
So my--
AUDIENCE: Spatial locality.
PROFESSOR: I'm sorry?
You don't have any spatial locality here.
AUDIENCE: I'm just wondering why they can't hold units.
I guess this is the question, why can't they hold multiple
addresses per cache line.
So why is it even pushed out?
It's being pushed out [UNINTELLIGIBLE] one per cache
line, right?
PROFESSOR: No, so it's getting pushed out because the cache
can hold M/B blocks, right?
So once it's accessed m over b blocks, if I want to access
anything else, something has to go out.
It's a capacity issue.
I access m over b blocks, something has to go out.
LRU says, the latest thing that I accessed, well that was
the first one, gets knocked out.
So what happens is every one causes a miss.
Even though I may access that very nearby in the future, it
doesn't take advantage of that.
Because LRU says knock it out.
AUDIENCE: [INAUDIBLE]
PROFESSOR: Is it row major that's the confusion?
This is the way we've been dealing with are matrices.
So in row major, there's a good--
that's nice, there's no chalk here.
Oh, there's a big one there, great.
Yeah, so here's B. So the order that B is stored is like
this in memory.
So basically we're storing these elements in this order.
So it's a linear block of memory, right?
And it's being stored row by row as we go through.
So actually if I do it like this, let me do this
a little bit more.
So the idea is that the first element is going to be here,
and then we get up to n minus 1, and then we get to n minus
2 is stored here.
n plus 1, n plus 2, 2n minus 1.
So that's the order that they're stored
in linear and memory.
Now these guys will all be on the same cache line if it's in
row-major storage.
So when I'm accessing B, I'm going and I'm accessing zero,
then I'm accessing the thing at location n.
And I'm going down like this.
At some point here I reach the limit of my cache.
This is M/B. Notice it's a different B--
script B verses--
so when I get to m over b.
Now all these things are sitting in
cache, that's great.
However, now I go to one more, and it says OK, all those
things are sitting in cache, which one do I kick out?
And the answer is the least recently used one.
That's this guy goes out.
AUDIENCE: Do you only use one element per cache link?
PROFESSOR: And I've only used one element from each cache
line at this point.
Then I go to the next one and it knocks out the second one.
By the time I get to the bottom and then I go up to the
top to access 1 here, it's not in cache.
And so it repeats the same process, missing
every single time.
We have a question.
AUDIENCE: Yeah,so my question is why does the cache know
where each row is?
To us, we draw the matrix, but the computer
doesn't know it's a matrix.
To the computer, its a linear array of numbers.
PROFESSOR: That's correct.
AUDIENCE: So why would it load the first couple elements in
the first row, and the second column is an extended row the
second time.
PROFESSOR: So the cache blocks are determined by the locality
and memory.
AUDIENCE: So my assumption would be the first
cache line for say--
PROFESSOR: Let's say 0 through 3.
AUDIENCE: So yeah, 0 through 3.
PROFESSOR: Let's say we have four items on that cache line.
AUDIENCE: 4 to 6.
PROFESSOR: The next one the hold 4 to 7, I think.
AUDIENCE: 4 to 7, yeah.
So that is a--
PROFESSOR: So that would be the next one, right?
4 to 7.
AUDIENCE: When you get cache line.
You are not using the fully cache line.
There is no spatial locale.
You are using one from the cache line.
PROFESSOR: So this code is using this one
and then this one.
It's not using the rest.
So it's not very efficient code.
AUDIENCE: So the cache line is holding the 0 to 3
and the 4 to 7.
[INAUDIBLE]
n plus 2 just reading--
[INTERPOSING VOICES]
PROFESSOR: Right.
And those are fixed.
So if you just a dice up memory in our machine into 64
byte sizes, those are the things that come in together
whenever you access anything on that line.
AUDIENCE: And on this particular axis, you never
actually get the 4 to the 7 in the--
PROFESSOR: Well we eventually do.
Until we get there, yes, that's right.
Until we get there.
Now, of course, we're also accessing A and C, but turns
out to this analysis it sufficient to show that we're
getting n cubed misses just on the matrix B. In order to say
hey, we've got a lot of misses here.
So this was the case where n was bigger than
the size of a cache.
So the situation is a little bit different if n is large
but still actually less than m over b.
So in this case, we suppose it n squared is bigger than m.
So the matrix doesn't fit in memory.
So that's what this part of the equation is, m to the 1/2
less than n is the same as n squared is bigger than memory.
So we still don't fit in memory, but in fact it's less
than some constant times m over b.
And now let's look at the difference with what happens
with the caches as we go through the algorithm.
So we essentially do the same thing.
Once again, we're going to assume LRU.
And so what happens is we're going to go down
a single row there.
But now, notice that by the time I get down to the bottom,
basically I've accessed fewer than some constant times m
over b locations.
And so nothing has gotten kicked out yet.
So when I go back to the top for the next access to B, all
these things are still in memory.
So I don't take a cache fall, a cache miss in those cases.
And so we keep going through.
And basically this is much better because we're actually
getting to take advantage of the spatial locality.
So this algorithm takes advantage
of the spatial locality.
If n is really big it doesn't, but if n is just kind of big,
then it does.
And then if n is small enough, of course, it all fits in
cache and there's no misses other than those needed to
bring it in once.
And then the same thing happens once you go through
the next one.
So in this case, what's happening is we have n squared
over b misses per run through the matrix B, and then we have
n times that we go through.
Once for every row of A. So the total then is n cubed over
b the cache block size.
So depending upon the size, we can analyze with this, that
this is better because we get a factor of B improvement.
But it's still not particularly good.
And it only works, of course, if my side of my matrix fits
in the number of lines of cache that I have.
Yeah, question.
AUDIENCE: Can you explain in--
I don't understand why you have n cubed over b?
PROFESSOR: OK, so we're going through this matrix n times.
And for each one of those, we're running
through this thing.
So this thing basically, I get to go b times through, because
all these things are going to be in memory when I come back
to do them again.
And so it's only once every B columns that I take a miss.
I take a miss and then I get to the other b minus 1 access
is that I get cache hit.
And so the total here is then n squared over b.
So therefore a total of n cubed over b.
So even this is not very good compared to what we can
actually do if we exploit the cache well.
So let's go on and take a look.
We saw this before.
Let's use tiling.
So the idea of tiling is to say, let's break our matrix
into blocks of s times s size.
And essentially what we do is we treat our big matrix as if
we're doing block matrix multiplications of things of
size s by s.
So the inner loop here is doing essentially the matrix
multiplication.
It's actually matrix multiply and add.
The inner three loops are just doing ordinary matrix
multiplication, but on s-sized matrices.
And the outer loop is jumping over matrix by matrix for each
of those doing a matrix multiply as
its elemental piece.
So this is the tiling solution that you've seen before.
We can analyze it in this model to see,
is this a good solution.
So everybody clear on what the code does?
So it's a lot of four loops, right?
Yeah.
AUDIENCE: There should be less than n somewhere?
There's like an and something.
PROFESSOR: Oh yeah.
That must have happened when I coped it.
That should be j less than n here.
It should just follow this pattern.
i less than n, k less than n, that should be
j less than n there.
Good catch.
I did execute this.
That must have happened when I was editing.
So here's the analysis of work.
So what's going on in the work?
So here we have, basically the outer loop is going n over s
times, each loop.
So there's cube there times the inner loops here which are
going each s times.
So times s cubed.
Multiply that through, n cubed operations.
That's kind of what you'd expect.
What about cache misses?
So the whole idea here is that s becomes a tuning parameter.
And whether we choose s well or poorly influences how well
this algorithm works.
So the idea here is we want tune s so that the submatrices
just fit into cache.
So in this case, if I want a matrix to fit into cache, I
want to be about the size of the square root
of the cache size.
And this is where we're going to use the tall-cache
assumption now.
Because I want to say, it fits in cache, therefore I can just
assume it all fits in cache.
It's not like the size fits but the actual data doesn't,
which is what happens with the short cache.
So the tall-cache assumption implies that when I'm
executing one of these inner loops, what's happening?
When I'm executing one of these linear loops, all of the
matrices are going to fit in cache.
So all I have is my cold misses, if
any, on that submatrix.
And how many cold misses can I have?
Well the size of the matrix is s squared and I get to bring
in b bytes of the matrix each time.
So I get s squared over b misses per submatrix.
So that was a little bit fast, but I just want to make sure--
it's at one level straightforward, and the other
level it's a little bit fast.
So the point is that the inner three loops I can analyze if I
know that s is fitting in cache.
The inner three loops I can analyze by saying, look it's s
squared data.
Once I get the data in cache, if I'm using an optimal
replacement, then it's going to stay in there.
And so it will cost me s squared over b misses to bring
that matrix in for each of the three matrices.
But once it's in there, I can keep going over and over it as
the algorithm does.
I don't get any cache misses.
Because those all fitting in the cache.
Question?
Everybody with me?
OK.
So then I basically have the outer three loops.
And here I don't make any assumptions whatsoever.
There's n over s iterations for each loop.
And there's three loops.
So that's n over s cubed.
And then the cost of the misses in the inner loop is s
squared over b.
And that gives me n cubed over bm to the 1/2 if you plug in s
being m to the 1/2.
So this is radically better because m is usually big.
Especially for a higher level cache, for an L2 or an L3. m
is really big.
What was the value we had before for the best case for
the other algorithm when it didn't fit in matrix?
It was n cubed over b.
b is like 64 bytes.
m is like the small L1 cache is 32 kilobytes.
So you get to square root the 32 kilobytes.
What's that?
So that's 32 kilobytes is 2 to the 15th.
So it's 2 to the 7.5.
2 to the 7 is 128.
So it's somewhere between 128 and 256.
So if we said 128, I've got a 64 and a 128 multiplier there.
Much, much better in terms of calculating.
In fact, this is such that if we tune this properly and then
we say, well what was the cost of the cache misses here,
you're not going to see the cost of the cache misses when
you do your performance analysis.
It's all going to be the work.
Because the work is still n cubed.
The work is still n cubed, but now the misses are so
infrequent, because we're only getting one every--
on the order of 64 times 128, which is 2 to the 6th times 2
to the 7th is 2 to the 13th is 8K.
Every 8,000 or so accesses there's a constant factor in
there or whatever, but every 8,000 or so accesses we're
getting a cache miss.
Uh, too bad.
If it's L1, that cost is four cycles rather than one.
Or that cost us 10 cycles if I had to go to L2 rather than
one, or whatever.
So is to the point is, that's a great multiplier to have.
So this is a really good algorithm.
And in fact, this is the optimal behavior you can get
for matrix multiplication.
Hong and Kung proved back in 1981 that this particular
strategy and this bound was the best you could do for
matrix multiplication.
So that's great.
I want you to remember this number because we're going to
come back to it.
So remember it's b times n to the 1/2, b times square root
of m in the denominator.
Now there's one hitch in this story.
And that is, what do I have to do for this
algorithm to work well?
It says right up there on the slide.
Tune s.
I've got to tune s.
How do I do that?
How do I tune s?
How would you suggest we tune s?
AUDIENCE: Just run a binary [INAUDIBLE].
PROFESSOR: Yeah, do binary search on s to find out what's
the best value for s.
Good strategy.
What if we guess wrong?
What happens if, say, we tune s, we get some value for it.
Let's say the value is 100.
So we've turned it.
We find 100 is our best value.
We run it on our workstation, and somebody else has another
job running.
What happens then?
That other job starts sharing part of that cache.
So the effective cache size is going to be smaller than what
we turned it for.
And what's going to happen?
What's going to happen in that case?
If I've tuned in for a given size and then I actually have
to run with something that's effectively a smaller cache,
does it matter or doesn't matter?
AUDIENCE: Is it still tall?
PROFESSOR: Still tall.
AUDIENCE: [INAUDIBLE]
PROFESSOR: So if you imagine this fit exactly into cache,
and now I only have half that amount.
Then the assumption that these three inner loops is running
with only s squared over b misses is going to be totally
out the window.
In fact, it's going to be just like the case of the first
algorithm, the naive algorithm that I gave.
Because the size of matrix that I'm feeding it, s by s,
isn't fitting in the cache.
And so rather than it being s squared over b accesses, it's
going to be much bigger.
I'm going to end up with essentially s cubed accesses
if the cache, in fact, gets enough smaller.
It's also one thing you have to put in there is what I like
to call voodoo.
Whenever you have a program and you've got some parameters
that, oh good we've got these parameters we get to tweak to
make it go better.
I call those voodoo parameters.
Because typically setting them is not straightforward.
Now there are different strategies.
One, as you say, is to do binary search by doing it.
There are some programs, in fact, which when you start
them up you call an initialization routine.
And what they will do is automatically check to see
what size is my cache and what's the best size should I
do something on and then use that when you actually run it
later in the program.
So it does an automatic adaptation
automatically when you start.
But the more parameters you get, the more
troublesome it becomes.
So let's take a look.
For example, suppose we have a two-level cache rather than a
one-level cache.
Now I need to have something that I tune in for L1 and
something that I tune for L2.
So it turns out that if I want to optimize s and t, I can't
do it in more with binary search because I have two
parameters.
And binary search won't suffice for figuring out
what's the best combination of s and t.
And generally multidimensional searches are much harder than
one-dimensional searches for optimizing.
Moreover, here's what the code looks like.
So now I've got, how many four loops?
1,2,3,4,5,6,7,8,9 nested for loops.
So you can see the voodoo is starting to
make this stuff run.
You really have to be a magician to tune these things
appropriately.
I mean, if you can can do it that's great.
But if you don't do it, OK.
So now what about three levels of cache?
So now we need three tuning parameters.
Here s, t and u, we have 12 nested four loops.
I didn't have the heart to actually write out the code
for the 12 nested for loops.
That just seemed overhead.
But our new halo machines, they have
three levels of caches.
So let's tune for all the levels of caches.
And as we mentioned, in a multi-program environment, you
don't actually know what the cache size is, what other
programs are running.
So it's really easy to mistune these parameters.
AUDIENCE: [INAUDIBLE] don't you have a problem because
you're running the program for a particular n, and you don't
necessarily know whether your program is going to run faster
or slower--
PROFESSOR: Well what you're usually doing, is you're
tuning for s not n, right?
So you're assuming--
AUDIENCE: No, no a particular n.
PROFESSOR: But the tuning of this is only dependent on s.
It doesn't depend on n.
So if you run it for a sufficiently large n, I think
it's reasonable to assume that the s you get would be a good
s for any large n.
Because the real question is, what's fitting in cache?
Yeah--
AUDIENCE: How long does it take to fill up the cache
relative to the context each time?
PROFESSOR: Generally you can do it pretty quickly.
AUDIENCE: Right.
So why does it matter if your have multiple users, if you
can fill it [INAUDIBLE].
PROFESSOR: No because, he may not be using all
of the cache, right?
So when you come back, you're going to have it polluted with
a certain amount of stuff.
I think it's a good question.
AUDIENCE: [INAUDIBLE]
PROFESSOR: OK, so anyway, so this is the-- yeah question.
AUDIENCE: So if n is really large, is it possible that the
second row of the matrix never loaded?
PROFESSOR: If n is really large--
AUDIENCE: Because n is really large, right?
Just the first row of the matrix will
fill up all the caches.
PROFESSOR: It's LRU and in B you're going down this way.
You're accessing things going down.
OK, good.
So let's look at a solution to these alternatives.
What I want to in particular take a look at is recursive
matrix multiplication.
So the idea is you can do divide and conquer on
multiplying matrices because if I divide each of these into
four pieces, then essentially I have 8 multiply adds of n
over 2 by n over 2 matrices.
Because I basically do these eight multiplies each going
into the correct result.
So multiply A11 B11 and add it into C11.
Multiply A12 B21 add it into C11 and so forth.
So I can basically do divide and conquer.
And then each of those I
recursively divide and conquer.
So what's the intuition by why this might a
good scheme to use?
AUDIENCE: [INAUDIBLE]
PROFESSOR: Well, we're not going to do parallel yet.
Just why is this going to use the cache well?
AUDIENCE: [INAUDIBLE]
PROFESSOR: Yeah eventually I get down to a size where the
matrix that I'm working on fits into cache, and then all
the rest of the operations I do are all
going to be cache hits.
It is taking something and it it's doing what the tiling is
doing but doing it blindly.
So let's take a look.
Here's the recursive code.
So here I have the base case if n is 1, I basically have a
one by one matrix and I just simply update c,
with a times b.
And otherwise what I do, is I'm going to do this by
computing offsets.
So generally when you're dealing with matrices,
especially if you want fast code, I usually don't rely on
two-dimensional addressing, but rather do the addressing
myself and rely on the compiler to do common
subexpression elimination.
So, for example, here what I'm going to do
is compute the offsets.
So here's how I do it.
So first of all, in practice what you do, is you don't go
down to n equals 1.
You have some cutoff.
Maybe n is 8 or something.
And at that point you go into a specialized routine that
does a really good 8 by 8 multiply.
And the reason for that is you don't want to have the
function call overheads.
This function call is expensive to do two floating
point operations here.
So you'd like to have a function call and then do 100
floating point operations or something.
So that you get a better balance.
Do people understand that?
So normally to write recursive codes you want a
course in the recursion.
Make it so you're not going all go the into way
down to n equals 1.
But rather are stopping short and then doing something that
doesn't involve a lot of overhead in the base case of
your recursion.
But here I'll explain it as if we went all the way
down to n equals 1.
So then what we do is, if this is a submatrix, which is
basically what I'm showing here.
We have an n by n submatrix.
And it's being pulled out on a matrix of size row size, of
width row size.
So what I can do is, if I want to know where the elements of
the beginning of matrices are, well the first one is exactly
the same place that the input matrix is.
The second one is basically I have to add n over 2 to the
location in the array.
The third one here, 21, I have to basically add n over 2 rows
to get the starting point of that matrix.
And for the last one I have to add n over 2 and n over 2 rows
and n over 2 plus 1 rows to get to that point.
So I compute those and now I can recursively multiply with
sizes of n over 2 and perform the program recursively.
Yeah--
AUDIENCE: So you said it rightly.
You're blindly dividing the matrix up til you get
something to fit the cache.
So essentially--
PROFESSOR: Well and you're continuing.
The algorithm is completely blind all the way
down to n equals 1.
AUDIENCE: This could never be better if the other one-- your
computer's version is well-tuned.
Because the applications are the same, but this one you
have all the overhead from the [INAUDIBLE].
PROFESSOR: Could be.
AUDIENCE: At the end, you still need to make a
multiplication and then go back and look at all of the--
PROFESSOR: Could be.
So let's discuss that later at the end when we talk about the
differences between the algorithms.
Let's at this point, just try to understand what's going on
in the algorithm.
Question--
AUDIENCE: [INAUDIBLE]
PROFESSOR: n over 2 times row size plus--
plus n over 2.
It should be row size plus 1.
You're right.
Good, bug.
Should be n over 2 times row size plus 1.
So let's analyze the work assuming the code
actually did work.
So the work we can write a recurrence for.
So here we have the work to solve an
n by n matrix problem.
Well if n is 1, then it's just order one work--
constant amount of work.
But if n is bigger than 1, then I'm solving eight
problems of size n over 2, plus doing a constant amount
of work to divide all those up.
So everybody understand where I get this recurrence?
Now normally, as you know, when you do algorithmic work,
we usually omit this first line because we assume a base
case of constant if it's one.
I'm actually going to keep it.
And the reason is because when we do caching the basic cases
are important.
So everybody understand where this recurrence came from?
So I can use the master theorem or something like that
to solve this.
In which case the answer for this is what?
Those of you who have the master
theorem in your hip pocket.
What's the solution of this recurrence?
People remember?
Who's has heard of the master theorem?
I thought that was kind of a prerequisite or something of
this class, right?
So you might want to brush up on the master theorem for the
quiz next week.
So basically it's a n over b, so it's n to the
log base 2 of 8.
So that's n cubed, n to the log base 2 of 8
is n to the n cubed.
And that's bigger than the order one here, so the answer
is order n cubed.
Which is a relief, right?
Because if weren't order n cubed we would be doing a lot
more work than one of the looping algorithms.
However, let's actually go through and understand where
that n cubed comes from.
And to do that I'm going to use the technique of a
recursive tree, which I think all of you have seen.
But let me go through it slowly here to make sure,
because we're going to do it again when we do cache misses
and it's going to be more complicated.
So here's the idea.
I write down the left hand side the recurrence, w of n.
And now what I do is I substitute, and I
draw it out as a tree.
I have eight problems of size n over 2.
So what I do is I replace that with the thing that's on the
right hand side, I've dropped the theta here, but basically
put just a constant one here.
Because I'll take into account the thetas at the end.
So I have a one here, and then I have, loops
that should be a w.
Should be w n over 2.
That's a bug there.
And then I replace each of those.
OK, wn over 2, sorry that should be wn over 4.
Ah, more bugs.
I'll fix them up on after lecture.
So this should be w of n over 4.
And we go all the way down to the bottom to where I hit the
base case of theta 1.
So I built out this big tree that represents, if you think
about it, that's exactly what the algorithm is going to do.
It's going to walk this tree doing the work.
And what I've just simply put up here is to work it does at
every level.
So the first thing we want to do is figure out what's the
height of this tree.
Can somebody tell me what the height of the tree is?
It is a log n.
What's the base?
Log base 2 of n, base because at every level if I hadn't
made a mistake here, I'm actually having the argument.
So I'm having the argument at each level.
So the height is log base 2 of n.
So LG is notation for log base 2.
So if I have log base 2 of n, I can count how many leaves
there are to this tree.
So how many leaves are there?
Well I'm branching a factor of eight at every level.
And if I'm going log base 2 levels, the number of leaves
is 8 to the log base 2.
So 8 to the log base 2 of n.
And then with a little bit of algebraic magic that turns out
that's the same as n to the log base 2 of 8.
And that is equal to n cubed.
So I end up with n cubed leaves.
Now let's add up all the work that's in here.
So what I do is I add across the rows.
So the top level I've got work of one.
The next level I work of eight.
The next I have work of 64.
Do people see the pattern?
The work is growing how?
Geometrically.
And at this level I know that if I add up all the leaves
I've got work of n cubed.
Because I've got n cubed leaves, each of
them taking a constant.
And so this is geometrically increasing, which means that
it's all born in the leaves.
So the total work is order n cubed.
And that's nice.
It's the same work is the looping versions.
Because we don't want to increase that.
Questions?
Because now we're going to do cache misses and it's going to
get hairy, not too hairy, but hairier.
So here we're going to cache misses.
So the first thing is coming up with a recurrence.
And this is probably the hardest part, except for the
other hard part which is solving the recurrence.
So here what we're doing is, we have the same thing is that
I'm solving eight problems of size n over 2 and to do the
work in here.
I'm taking basically order one cache misses.
However I do, those things work out.
Plus the cache misses I have in there.
But then at some point, when I'm claiming is that I'm going
to bottom out the recursion early.
Not when I get to n equals 1, but in fact when n squared is
less than some constant times the cache size.
For some sufficiently small concept.
And what I claim, at that point, is that the number of
cache misses I'm going to take at that point, I can just,
without doing any more recursive stuff, I can just
say it's n squared over b.
So where does that come from?
So this basically comes from the tall-cache assumption.
So the idea is that when n squared is less than a
constant times the size of your cache, constant times the
size of m, then that means that this fits into--
the n by n matrices fit within m.
I've got three of them.
I've got C, A and B. So that's where I need a constant here.
So they're all going to fit in memory.
And so if I look at it, all I have to do is count up the
cold misses for bringing in those submatrices at the time
that n hits this threshold here of some constant times m.
And to bring in those matrices is only going to cost me n
squared over b cache misses.
And once I've done that, all of the rest of the recursion
that's going on down below is all operating out of cache.
It's not taking any misses if I have an
optimal replacement algorithm.
it's not taking any more misses as I get further down.
Questions about this part of the recurrence here?
So people with me?
So when I get down to something of size n squared,
where the submatrix is size n squared, the point is that
I'll bring in the entire submatrix.
But all the stuff that I have to do in there is never going
to get kicked out, because it's small
enough that it all fits.
And an optimal algorithm for replacement is going to make
sure that stuff stays in there, because there's plenty
of room in the cache at that point.
There's room for three matrices in the cache and a
couple of other variables that I might need and that's
basically it.
Any questions about that?
So let's then solve this recurrence.
So we're going to go about it very much the same way.
We make draw a recursion tree.
So those of you are rusty in drawing recursion trees, I can
promise you there will be a recursion tree on the quiz
next Thursday.
I think I can promise that.
Can I promise that?
Yeah, OK I can promise that.
The way I like to do it, by the way, is not to try to just
brought out all at once.
In my own notes when I do this I always draw it step by step.
I copy over and just do a step by step.
You might think that that's extensive.
Gee, why do I have to draw every one along the way?
Well the answer is, it's a geometric process.
All the ones going up to the last one are a small amount of
the work to draw out the last one.
And they help you get it correct the first time.
So let me encourage you to draw out the tree
iteration by iteration.
Here I'm going to just do replacement.
So what we do is we replace with the right hand side to do
the recursion.
And replace that.
And once again I made the bug, that should be in over 8.
Sorry, n over 4 here.
n over 4.
And then we keep going down until I get to the base case,
which is this case here.
Now comes the first hard part.
How tall is this tree?
Yeah--
AUDIENCE: [INAUDIBLE]
square root of n over b.
You want n squared to be cm, not [INAUDIBLE].
PROFESSOR: So here's the thing, let's discuss, first of
all, why this is what it is.
So at the point where n squared is less than cm, that
says that it's going to cost us n squared over b.
But n squared is just less than cm, so therefore, this is
effectively m over b.
Good question.
So everybody see that?
So when I get down to the bottom, it's basically costing
me something that's about the number of lines I have in my
cache, number of misses to fill things up.
The tricky thing is, what's the height?
Because this is crucial to getting this kind of
calculation right.
So what is the height of this tree?
So I'm having every time.
So one way to think about it is, it's going to be log bas 2
of n, just as before, minus the height of the tree that is
hidden here that I didn't have to actually go into because
there are no cache misses in it.
So that's going to occur when n is approximately m, cm,
sorry when n is approximately square root of cm.
So I end up with log of n minus 1/2 log of cm.
That's the height here.
Because the height at this point of the tree that's
missing because they're no cache, I don't have to account
for any cache misses in there, is log of cm to the one half,
based on this.
Does that follow for everybody?
People comfortable?
Yeah?
OK, good.
So now what do we do?
We count up how many leaves there are.
So the number of leaves is 8, because I have a branching
factor of 8, 2 whatever the height is.
Log n minus 1/2 log of cm.
And then if I do my matrix magic, well that part is n
cubed, the minus becomes a divide, and now 8 to the 1/2
log of cm is the square root of n cubed,
which is m to the 3/2.
Is that good?
The rest of it is very similar to what we did before.
At every level I have a certain number of things that
I'm adding up.
And on the bottom level, I take the cost here, m over b,
and I multiply it by a number of leaves.
When I do that I get, what?
I get n cubed over b times m to the 1/2.
This is geometric.
So the answer is going, in this case, just going to be
the sum of a constant factor times the large thing.
And why does this look familiar?
That was the optimal result we got from tiling.
But where's the tuning parameters?
No tuning parameters.
So that means that this analysis that I did for one
level of caching, it applies even if you have three levels
of caching.
At every level you're getting near optimal cache behavior.
So it's got the same cache misses as with tiling.
These are called cache-oblivious algorithms.
Because the algorithm itself has no tuning parameters
related to cache.
Unlike the tiling algorithm.
That's a cache-aware algorithm.
The cache-oblivious algorithm has no tuning parameters.
And if it's an efficient one.
So, by the way, our first algorithm was
cache-oblivious as well.
The naive one.
It's just not efficient.
So in this case we have an efficient one.
It's got no voodoo turning of parameters, no explicit
knowledge of caches, and it passively autotunes itself.
As it goes down when it fits things into cache it fits them
and uses things locally.
And then it goes down and it fits into the next level of
cache and uses things locally and so forth.
It handles multi-level caches automatically.
And it's good in multi-programmed environments.
Because if you end up taking away some of the cache it
doesn't matter.
It still will end up using whatever cache is available
nearly as well as any other program could use that cache.
So these are very good in multi-programmed environments.
The best cache-oblivious matrix multiplication, in fact
doesn't do an eight way split as I described here.
That was easier to analyze and so forth.
The best one that I know work on
arbitrary rectangular matrix.
And what they do, is they do binary splitting.
So you would take your matrix, i times j, So if you take a
matrix, let's say it's something like this.
So here we have i, k, k, j.
And you're going to get something of shape.
i times j, right?
What it does, is it takes whatever
is the largest dimension.
In this case k is the largest dimension.
And it partitions either one or both of the
matrices along k.
In this case, it doesn't do that.
And then it recursively solves the two
sub-rectangular problems.
And that ends up being a very, very efficient fast code if
you code that up tightly.
So it does binary splitting rather than--
and it's general.
And if you analyze this, it's got the same behavior as the
eight way division.
It's just more efficient.
So questions?
We had a question about now comparing
with the tiled algorithm.
Do you want to reprise your question?
AUDIENCE: What I was saying was, I guess
this answers my question.
If you were to tune the previous algorithm properly,
and you're assuming it's not in a multi-program
environment, the recursive one, it will never be the one
that is locked.
[INAUDIBLE]
PROFESSOR: So at some level that's true, and at some level
it's not true.
So it is true in that if it's cache-oblivious you can't take
advantage of all the corner cases that you would might be
able to take advantage of in a tiling algorithm.
So from that point of view, that's true.
On the other hand, these algorithms work even as you go
into paging and disks and so forth.
And the interesting thing about a disk, if you start
having a big problem that doesn't fit in memory and, in
fact, is out of core as they call it, and is paging to
disk, is that the disk sizes of sectors that can be brought
efficiently off of a disk, vary.
And the reason is because in a disk, if you read a track
around the outside you can get two or three times as much
data off the disk as a track that you read near the inside.
So the head moves in and out of the disk like this.
It's typically on a pivot and pivots in and out.
If it's reading towards the inside, you get blocks that
are small versus blocks that are large.
This is effectively a cache line size that
gets brought in.
And so the thing is that there are actually programs in
which, when you run them on disk, there is no fixed size
tuning parameter that beats the cache-oblivious one.
So the cache-oblivious one will beat every fixed-size
tuning parameters you put in.
Because you don't have any control over where your file
got laid out on disk and how much it's bringing in and how
much it isn't varies.
On the other hand, for in-core thing, you're exactly right.
That, in principle, you could tune it up more if you make it
more cache aware.
But then, of course, you suffer from portability loss
and from, if you're in a multi-programmed environment
and so forth.
So the answer is, that there are situations where you're
doing some kind of embedded or dedicated type of application,
you can take advantage of a lot of things that you want.
There are other times where you're doing a
multi-programmed environment, or where you want to be able
to move something from one platform to another without
having to re-engineer all of the tuning and testing.
In which case it's better to use the cache oblivious.
So as I mentioned, my view of these things is that
performance is like a currency.
It's a universal medium of exchange.
So one place you might want to pay a little bit of
performance is to make it so it's very portable as for the
cache-oblivious stuff.
So you get nearly good performance, but now I don't
have that headache to worry about.
And then sometimes, in fact, it actually does as well or
better than the one.
For matrix multiplication, the best algorithms are the cache
oblivious ones that I'm aware of.
AUDIENCE: [INAUDIBLE]
currency and all the different currencies.
Single currency.
PROFESSOR: You want a currency for--
so in fact the performance for this is people who have
engineered it to take advantage of exactly the cache
size, we can do just as well with the cache oblivious one.
And particularly, if you think about it, when you've got
three levels hierarchy, you've got 12 loops.
And now you're going to tune that.
It's hard to get it all right.
So next time we're going to see a bunch of other examples
of cache-oblivious algorithms that are optimal in terms of
their use of cache.
Of course, by the way, those people who are familiar with
Strassen's algorithm, that's a cache-oblivious algorithm.
Takes advantage the same kind of thing.
And in fact you can analyze it and come up with good bounds
on performance for Strassen's algorithm just the same.
Just as we've done here.