the second section of course of course wason a applications in geometry in
approximations statistics areas like that. And this last section in the course is on.
Numerics and algorithms how do youactually solve these problems?
Now, you have been solving in theseproblems more well for the
last two parts of the course you have beensolving these problems.
but presumably
at this point you don't know how it's doneand the main
goal of this section, the whole section ofthe course is to demystify.
How it is that we solve these convexproblems, right.
So that's that's the idea.
Now, most of you are not going to writeyour own solvers.
Right.
Actually, some of you might have to, ifyou solve, if
you go into areas where the problems areabsolutely gigantic, you know,
something like CVX is not going to work,
you're going to have to develop your ownsolver.
or if you go into like real time systemswhere you
need to be able to solve optimizationproblems in micro seconds.
Again something like cvx is not going towork for you
and you're going to really have to knownhow its done.
I should say a couple more things aboutthat.
I still believe strongly, its not just aquestion of demystifying it
but I believe its very important foreveryone to know a little bit
About how these problems are solved andspecifically what you want to know
is the link between problem structure andhow fast we can solve it.
So that's really, if there was one takeaway for
the entire last part of this course, it'sto understand that.
That one thing.
You should be able to look at that problemand say "Ooh!
That's got a hundred thousand variables, but I can solve that super fast, because
I recognize the structure." So that's whatthis whole section of the course is about.
Now we are going to go back today andstart, well, at the very beginning.
We are going to talk a little aboutnumerical algebra.
Some of you probably know this, others ofyou probably don't.
People who have used things like...
Matlab have certainly been using this as auser, as
a consumer but you may not know actuallyhow it's done.
and what we'll do today is just look atthe basics
and I should also say how it fits into thebig picture.
So the big picture is, this lecture, we'regoing to
talk about How do you solve a linearequation.
So that's fine.
And there's a finite algorithm for thatand we'll look at
how that's done, and we'll we'll see we'llsee how it's done.
It's gaussian elimination.
You've been doing this since high school Iassume.
Okay then I'll tell you how this fits intothe rest of the course.
Right?
So in the rest of the course we'll see howdo you
minimize, we'll start by talking how doyou minimize a smooth function.
well, with start with how do you, we'll
start with how do you minimize a quadraticfunction.
Well, you set the gradient equal to 0.
But the grading of a quadratic fuction isa linear
function, and so, minimizing a quadraticfunction, solving linear equations.
So if you like, you can think of thiswhole stuff we're going to look at now,
how do you solve linear equations, as howdo you minimize a quadratic funtion.
So, that's the optimization connection. Okay so, we will first start
by talking matrix structure and algorithmcomplexity.
matrix structure is, you know, somethingsare obvious like
what is, you know, diagonal matrices, identity matrices sparse matrices.
We'll talk about that.
then we'll talk about the dominant methodused to solve linear equations.
and that's using factored matrices.
So you'll factor a matrix and then solve,do a sequence
of solves, with the factored, with thefactors.
we'll talk about some very famousfactorizations, the so called
LU, that's going to be for lowertriangular, upper triangular, the
Cholesky factorization and somethingcalled
the LDL transpose factorization where L
is lower triangular, D is diagonal and Lis upper triangular.
Actually we will see in some cases D isblock diagram.
Then, we'll look at a general constructioncalled block elimination
in the matice inversion lemma and this issomething that shows you how
it is that you can solve
sometimes absolutely enormous sets oflinear equations.
very efficeintly by explointing somestructure that's present.
And the last thing we'll talk about,
though very briefly is solvingundetermined equations.
So first let's just start with solving Axequals b.
where A is n by n.
Okay?
And in fact let's say on that laptop rightthere,
that's not mine, but I know the, I knowthe number.
I wa, I want to ask you all, how long do
you think it takes to solve the equationAx equals b?
Where x is size 1,000.
Before you say anything, I want to point afew things out.
That's solving
1,000 equations, for 1,000 variables. The size of A, is a million coefficients.
Everybody got it?Okay?
And now what I want to know is how fastcan I solve it?
The answer is 30 milliseconds.
Okay, and I happen to know the number onthat thing is about ten.
So actually, this is an amazing thing, andit's actually
worthwhile for all of us to sit and thinkhow unbelievably amazing that is, right.
i mean that is ridiculous.
Right?
That's something that takes like abillion, we'll
see it takes, generic methods take abillion operations.
It's ten, fifteen miliseconds on that.Okay?
It's got multiple cores, does all sort ofthings are happening to make that number.
But just the idea actually that it's lessthan a second is unbelievably impressive.
Okay?So well we'll get to that.
So it, so it basically for general methodsright?
The general methods means you, you'reexploiting
nothing in the matrix a, absolutely none,right?
All the entries are assumed to be non-zerothere's no particular structure.
The structure would be something like amight be
A discreet foray transform or a discreetcosine transform.
In which case, of course, you can solve itfaster than n cubed.
You can solve it in like n log n. Okay so that would be structure.
Other structure would be sparsity.
Right, that would be that A has a lot ofzeroes.
Right, one definition of sparsity goessomething like this: it's.
A matrix that has enough zeros that it'sworth your while to take advantage of it.
That's 1967; I mean it's a perfectly gooddefinition, right?
So,
okay.
So general methods, it says that thesethings grow, like N cubed.
oh, so for fun, how long would it takethat laptop there to solve?
Say a hundred, hundred equations,ah and ahundred
variables that's 10 thousand coefficientsin the matrix.
So the answer was what 30? >> Micro seconds
>> Micros seconds, thank you, okay. So I mean that's ridiculous,
by the way this also hints Why you shouldknow about all those things, right.
These numbers are so ridiculous that itmeans these things can now be embedded
on all, it means the amount of computationyou can do can now, is unbelievable.
you can put stuff, you can actually embed, all of these things can be embedded.
They can, and things that would normallytake like seconds to run.
>> can run in milliseconds. They can be embedded in things.
I mean, crazy stuff
can happen.So, OK.
all right.
So, basically for general, these are thegeneric methods.
This grows as N cubed.
And I'll, I mean I'll tell you what thegeneric methods are, shortly.
what, what they are.
and this is, this is the part that Ithink.
I want everybody to know.
when A is structured, and we'll talk aboutall sorts of structure, it's way faster.
Or generally, substantially faster than ncubed.
Okay, and structured is stuff like banded.
A banded matrix is 1 Where, basically, a ij is
zero, if the difference between i and j isbig enough.
Bigger than some thresh hold.
Right, systems like that can actually besolved in linear time, okay.
So, I mean these are just things you did, you wouldn't I mean unless you knew this,
it's not obvious right, so what that sayis I can solve a million variables with
a million equations on that laptop, actually I can do that on my phone
but on that,that laptop I solve a
system of a million equations, a millionvariables.
If, if it's banded with a bandwidth oflike five.
Right?
And not only that, it can solve, I mean,unbelievably fast.
I won't do the calculation now but it'llbe shocking.
OK?
And if you didn't know that, you'd think
well that's, that's a big hard problem orsomething.
You multiply it by n cubed.Okay.
so, what are the basic ideas? Well this goes back to the 60s.
It's still sort of relevant.
It's less relevant that it was. I'll say a little bit about that.
But, it's still not a bad indicator.
so, basic concept is a FLOP, that'ssupposed to be a Floating Point Operation.
And roughly, it's something like anaddition, a
subtraction, a multiplication, division of2 floating point numbers.
Now, there's some older definitions where.
It's an addition plus a multiply and you
might include things like square roots inthere.
All of that is nonsense because squareroots cost like 50 times more.
A division is something like, cost you 20times more than an add.
So for all of these things, the, the keyidea now
is just to not to take the numbers tooseriously, right?
They're to be considered to be within afactor of ten.
Right?
And that, that, that's the right number.these,
the flop was invented sort of when, youknow,
main CPUs didn't do floating pointarithmetic, you know.
Then basically had to like send a messageover to
some giant board that was like yourfloating point board.
And all sorts of stuff happen.
And then when that multiplying wasfinished.
Right?Something would come back.
Everybody, you know, and so in thatsetting,
sure, counting the number of flops madetotal sense.
It's ridiculous now right because, youknow it, with the number
of, I mean there's all sorts of things youcould do there not even getting into
like minor details but you can do crazystuff all in parallel and stuff like that.
So, okay.
>> All right.
So, the, the, the idea, and this traces
to the 60s and is still actually somewhatrelevant.
it goes something like this.
To estimate the complexity of an
[UNKNOWN]
what you do, is you'd express the numberFLOPs.
it usually, it's a polynomial.
It's a polynomial function of the problemdimensions.
And then what you do is you drop. All but the leading terms, right?
So, you know, if it's something like, youknow, 1/2 n cubed plus 2 n
squared plus blah, blah, blah - somethinglike that - you just say, well, it's one.
The complexity's 1/3 n cubed.
And honestly, these days, the 1/3 doesn'tmake any difference at all.
Because you're not shooting for - youknow, factor of
2 is ridiculous.
You could not predict using these methods,right?
Factor of 10 is the right, the rightamount.
Now you know this is not an
accurate predictor of computation onmodern computers.
Right.
It, you can get it within a factor of ten
if, if thing you know in the absence ofcertain bad
things and actually if things are donereally well right
like it's optimized code it once againbecomes a reasonable predictor.
So it's kind of an interesting thing.
but, you know it's still, it's, it's veryuseful as a rough estimate of complexity.
It's still, it's still gives you that,right?
but, you know, do not come along and sayoh, you know, this
one is half the number of flops therefore, it must run faster or something.
I mean things like your, your accesspattern
and, and, and Locality of reference andmemory
are way more important than factors of to
and flop counts, and things like that, so,okay.
I mean another nice thing is that a lot ofthese things are
done for us by people who know whatthey're doing, namely those who
write like compilage and things like that, so that, that's the good news.
Okay, so.
let's look at some, operations and this isa conceptual division.
It's a vector-vector, matrix-vector, matrix-matrix product.
And I should say this is something youshould probably all know about.
there's no reason not to.
There's something called BLAS and ittraces to the 60s, 70s.
It's called it, it stands for Basic LinearAlgebra Subroutines.
Right?
And believe it or not, the, the standardsof these are written in I, well four tran.
That gives you a rough idea.
OK?And by the way still.
So, and then they have of course, theseare
re implimented in all sorts of reallanguages actually.
So anyway this is basic linear and theycall these
things, this is called blass- Level one,this is blast two,
and that's blast three.Right.
You don't need to know that, but actuallythat's not true.
An educated person just knows thesethings.
You don't have to know the details ofthese things.
But you should know but there's blastlevel
one, blast level two, and blast levelthree.
And I'll say a little bit about what, what
they mean and what why it's important toknow about.
Okay.
So vector-vector operations this would besomething
like, you know, calculating an innerproduct, right?
You calculate an inner product, yeah, youcan count the number of flops.
It's like 2n minus 1 flops, 2n if n islarge, right?
So that's something like that.
So you'd say it's order 2n or somethinglike that.
Another one would be like a sum.
Something, something like that.
by the way if you're taking an innerproduct it,
it what, would also depend if, if thevectors are sparse.
Then, needless to say this could be wayfaster.
And it actually depends, I mean,
quite a bit on You know, the datastructure
used to describe sparse vectors and thingslike that.
So.
Okay.
the other thing let's just say about aninner product, I mean just to give a
hint, sort of, about modern stuff, issomething
like an inner product, is, is easilyparalyzable.
Right so so on a modern on on some modernplatform it'd be
divided up if you do an inter plot of 2million long vectors
they have multiple threads and each onewould be taking inner product in
chunks and then those would be addedtogether in things like that so.
That's.
I'm just saying.
That's by the way how you get the 20millisecond,
30 millisecond numbers to solve a thousandequations on that laptop.
Okay?So, but it's not bad.
If you want to think of something seriallystepping
through like, please get me the next x i.
Please get me the next y i. Please multiply them together.
And then, plus equals that on some, some,you
know, you're inner product plus equals x iy i.
If you want to think of that. Go ahead, it doesn't hurt.
It's charming, it's simple.
But, you should understand that's not howit works.
I mean, unless you write the code rightlike that, right?
In which, in which case it'll run, whatever, 10 times slower, like so.
Okay.All right.
Next one is matrix vector product.
That's y equals a x.
so here you multiply a matrix by a vectorand what is it?
Well, you know there's many ways to thinkof it.
But, one way if you're just countingFLOPs, it doesn't matter, is that
it's basically, you take the inner productof each row of A with x.
Okay?
So, it's basically calculating M innerproducts in our N.
And that's something like 2 MM.
Right?
And you forget the 2. This is like,
that's a joke.Right?
It's order M, N.Right.
So, it's the number, it's vaguely numberentries in A.
now if A is sparse, this is obviously awhole lot less.
Right?And but actually,
[INAUDIBLE]
how fast it is, is quite a complicatedtopic and depends on exactly what
data structure's used to represent thesparse
matrix and all sorts of stuff like that.
Right?
But Again, it's not easy at least at a,without looking,
getting into the details to realize that,if a matrix is a
million by a million, with maybe 10million non-zero entries, you
only have to multiply it by the 10 millionnon-zero entries, right?
So, okay.and it can be less,
it's a structure.so, for example, here's one.
if I gave you the matrix A.
As a, a low rank product, right? So if I gave you the matrix A.
Let's suppose a is million by million. But I give it to you as a million by ten.
Multiply it by ten by million, right?
So it's what you would get if you ran pca,or sbd, or something
like the low rank approximation.
>> Number 1, you'd actually be able tostore it.
You can't store a million by millionmatrix, right?
So, but you could store it this way.
And the multiply would certainly not be amillion squared.
It would be way less.
It's order of a million.
In fact, it's like 20 million times thatbecause, 10 by 10.
Everybody got this?
Right.So, okay any of BLAS level three.
These are matrix matrix products.
So, that means multiply two matrices.
I mean that's the canonical one, ismultiply two matrices.
And here, it's actually the product of thedimensions, basically, right?
If these are like m by n, m by p.
It'd be, it's mnp.
So it's the product of the, of thedimensions.
Okay.
that are sparse, it's much less and youget different things
if if they're symmetric, you can dividethings by two, I mean,
that typically doesn't matter but that's,that's the idea.
and now, this is probably where I, I mean,
again, we, these are now things We're notgoing
into these details, but I'll just say alittle
bit about it so that you appreciate it,right.
So, everybody here knows how to multiplytwo matrices, right.
I mean it is, you could write a five-line
C program that would do this, every herecould, right.
It can be truly trivial. It could be like how do you get C i k?
Why, you, you walk along A and you walkalong
B and multiply the, the associatedentries.
And then add it all up.
right?
So, you could write a five line C program.Right?
It'd be a bunch of fours and wouldn't evenbe five lines, probably.
It'd be 5 lines with a lot of commentingbut, okay and you can compile it.
And, of course, it would do, like you knowwhatever.
Two MNP FLOPs by the time it runs.in fact.
if you were to do, if you were to use
LAPACK, that, that's the, again, this isjust for cultural background.
You do not need to know this although it'snot bad to know it, okay?
So, LAPACK is just, this is an open sourcepackage of linear algebra routines.
and, you know, multiplying two matrices.
This I mean doesn't seem like there be anyif you don't think about
this doesn't seem to be any subtlety thereright we all know what CIK
is, its sum over J, AIJ be JK right I meanit
don't look like a whole lot of room foroptimizing there right.
I can give you something a little
more complicated like sorting or somethinglike that
or finding the shortest path But then you
can imagine somebody smart could dosomething better.
Everybody, everybody see what I'm sayinghere?
I mean, come on, calculating the productof two matrices?
You know, one person's implementation, youcan beat a naive one by a factor of 100?
That doesn't seem right.
Does anyone know, by the way, how that'sdone, and why?
Its actually kind of interesting.
I mean, this is just for culturalbackground, but its good to know.
Here's how you multiply two thousand bythousand matrices.
Here's how you do it if you use thispackage.
What happens is, thousand by thousandmatrices
are broken down into things like ten
by ten block matrices with blocks that arehundred by hundred.
Everybody following this?
Now, matrix multiplication you can doblock wise.
Right.
So, what that says is you end upmultiplying 200 by 100 matrices
by blocking them and then each block youdo 100 by 100 multiply.
Okay.
Now the 100 by, how do you do 100 by 100multiply?
It's blocked down to 10.
So each of those is a ten by ten matrix often by ten matrices.
Everybody got this?
And then so what's really happening isyou're
doing giant piles of ten by tenmultiplicates, and
by the way the parameters of ten I made upare totally made up and it would differ.
Depends on your architecture.Alright?
Okay, so that's called blocking.
And the idea behind blocking is you canarrange it if you do it just right.
If the block's in the right sizes likethat ten by ten
multiplied might hap, might happen in, youknow, the ten by ten multiplied
might happen with registers or something,right.
And then the hundred by hundred might fitin L, might just fit in L1, right?
And stuff like that.
So if you do this right, if you do yourblocking right.
>> Then, you, and not only that.
You could do part of it in parallel andstuff like that.
If you have multiple cores, you know, somebody's multiplying these.
Somebody's multiplying these.
You exploit all this stuff.
And that's actually how matrices getmultiplied.
OK?
On,
at, at least on a, you know, on sort of amodern processor.
Right?
So, that, that, that's the idea.
And, and if you wa, and the memory accesspattern is such that there are never
[UNKNOWN].
And by the way, how do they, how do theychoose the block sizes?
That's they have something called atlas,right?
That's
[LAUGH]
really kind of cool.
I mean, you don't need to know any ofthese things, right?
But it's sort of fun.
That's automatically tuned linear algebrasoftware.
this all open source stuff.
And what it does, it takes your computerand it just said, I try 10, right?
And it does it, and it factors matricesfor a while
and it's say I'll try blocking it as 9 by9.
It tries 9, it tries 11, it tries 12, andit just does
bisection or just tries a bunch of themand then half an hour later.
It's decided on your architecture what the
right cash, what the right blocking sizesare.
Okay, completely imperdible, everybody gotthis?
So, and when you see these curves theirabsolutely amazing.
It would be things like, you would see
something going like, here's, here's afunction of n.
And say it would be, it would be lookinglike that.
That's a log, log plot.
Then all of a sudden you would see this.
[LAUGH]
Right, and, and you'd say, and by the way,what would this be?
Yeah, so this is something that's justbarely fitting in some cache, right?
And this is something that, that.
This is called bad blocking.
This is, this is where your block sizeswere exactly You
know a handful of things to big to fit insome cache.
So basically it's a disaster, somethinglike that.
In fact this can even go down that's thecool part.
This can be nonmonotonic but any way okay.
So that's more than you need to know aboutall of this but it is good to know.
actually it is good to know even
if you're interests are quite theoriticalor mathematical.
This is good to know beacuse if you dotheorical or very
mathematical things You know, here's, thisis how it hits the road.
It hits the road right here, in, in LAPack in most cases.
Well, except in giant sparse things, butwe'll talk about that.
But a lot of times, you do some coolstatistics, you do any you do
[UNKNOWN],
it's this.
Now we're going to talk about solvinglinear
equations, and The one way to do it.
We'll first talk about classes.
classes, of, of, matrices. Where solving Ax equals b is easy.
And then there's a simple trick that putsit all together.
So, one would be this.
diagonal matrices.
Right, so if A is diagonal then solving aequals b.
Is simple.
oh, by the way, one interesting thing is
people will, when people write A inverseb,
mathematically it's what you mean.That's what x is.
but it is universally understood, even ifpeople write it in an
algorithm statement that it doesn't, thatyou do not parse it this way.
A inverse times b because that wouldsuggest
code like, you know, a. in, you know, times
b and eh, eh, except in extremely rarecases you would never compute it that way.
But,.so, a inverse be is sort of like a, a
phrase that when someone describes analgorithm you should see.
And it doesn't mean form inverse and
multiply by b Mathematically it's thesame.
The semantics is identical.
But you should understand that this is socommon so, and people will even
say things like inverting that matrix is
easy and they don't mean inverting thematrix.
They mean solving the equation. Okay, so this is just, so, okay.
so here,
you know, computing a inverse b. If, if A is diagonal.
It's incredibly easy. It's just you divide.
It's N FLOPs. You know, roughly.
Whatever it is.OK?
the next one is interesting.Right?
It says that a matrix is lower triangular.Right?
And so, that says, you know it looks likethis.
A11, A21, A22, you know, something likethat.
Okay?
And you're multiplying this by, you knowx1
up to xn and that's supposed to be, youknow b1 up to bn.
Okay?
And yeah, how do you, how do you do this? You look at the first row.
The first row says a1, 1, x1 equals b1, sox1 is, is obviously a b1 array 1, 1.
Right?
Oh, and by the way, what exception do youthrow if a1, 1 is zero?
What if we have a divide by 0 exception?
Yeah.
You throw a message back which is like Asingular.
So, so I'm not going to calculate a Ainverse b for you.
Right, so okay.
So once you know x1, right? And there's lots of ways to think of this.
Once you know x1, the second row says A21x1 plus A22 x2 equals b two right?
But you know x one so you just plug inthat value subtract it on
the right and then divide it by A22.
So the algorithim looks osmething likethis an
dyou can sort of calculate the number offlops.
I mean it doesn't matter but you cankind of get a rough idea what it is.
>> This thing is kind of growing.
the, oh this, and it's got a, a greatname.
This is forward substitution.
And the reason it's called forwardsubstitution is you calculate XI.
And then XI, from then on, you substitute
that value into the equation.Right?
And so, you can see, you do kind of like.
Length 1, length 2, length 3, length 4,and it's growing.
And you know that the sum of those
numbers is something like the square over2, right?
And if you do it exactly right, it's justlike n squared, right?
You count blocks, that doesn't matter. It's n squared blocks.
Okay, and that's called forwardsubstitution.
And there's a, a similar one, if a matrixis upper-triangular.
And what do you imagine it would be, youknow, the same thing.
It's like called, it's backwardssubstitution.
Right?
You start by getting XN first and goingbackwards.
Okay?
So these are, and these are, justcompletely simple.
There's nothing complicated here.
this is really interesting when A issparse.
Because when you have a lower triangularsparse matrix.
I don't know.Go ahead and give me a rough number.
so for a sparse matrix.
1, and a very important attribute is
number of non zeros, sometimes writtenNNZ, okay?
So I just roughly, if I have, in terms ofthe number of non zeroes.
How many FLOPS
does it take to do forward or backward
substitution with a lower triangular, upper triangular matrix?
That's sparse, in terms of the number of
[INAUDIBLE]
zeros.What is it?
Yeah, it's winning.
In fact it equals, kind of, to the n andz.
Everybody got that?
Because, basically each of thesecoefficiants appears once.
Right?And so, it's about linear.
Right?
So, that means if I have a million bymillion matrice Right?
With let's say 50 million non zeros,that's you know
on average of 50 per non zeros per row orcolumn.
That's 50 million, 50 million flops.
How fast can I do it on there? This rough order magnitude.
You can assume that's between 1 and 10gig, gig gigaflops something like that.
So the answer is milliseconds.
another one would be orthogonal.
So if a matrix is orthogonal then in factsolving, well
computing a inverse b is the same as wella transpose b.
Right?
So its just a matrix vector multiplied, its ordered n squared in general
but it turns out it actually depends howthe orthogonal matrix is stored.
There are other data structures.
Here's a very famous one. And it looks like this.
That's, that's called, that's areflection.
Right?That's, that's an orthoganol matrix.
And you can see it is an identity minus arank one matrix.
Right?
We're going to talk a lot about that latertoday.
but so there that can be done in order nflops.
Right, so very fast.Right.
So, so for example, one, I mean, it's notcompletely standard, but
it, you would occur often is one datastructure for representing an orthogonal
matrix is as a product these things, thesereflections here.
And you would actually give them, youwould give the u , u1, u2, u3,
u4 and you can see the multiple, multipleactually would be the same to multiply,
multiply by the transpose which is the
inverse, each of those would simply belinear
in that in that data structure, right andif I now want this permutation matrix so
permutation matrix.Well, it's kind of silly.
Right?
If I solve, you know, P x equals b, all Ido is permute the interest, right?
And actually by this 1960's version ofthings.
When you say, how many flops is that, theanswer is zero.
There's no floating point in just doing,you're copying.
You're copying I, I guess you, all you'redoing is you're copying,
you know, whatever it is eight bytes or, how many bytes is a, a double?
Eight?
So you're copying
[UNKNOWN]
eight bytes from here to there and aroundlike that.
So, and oh, by the way, if that's ahundred
million long vector, and their accesspattern is completely random or
whatever This could, you could end upwith, this could be
not, it could be 0 flops and actually takesome time.
Okay?
[COUGH]
Alright okay, so alright and now
we get to, how do you solve Ax equals b? and, it's a general method.
It's, it's, it's, it's a high level, It's,it's, it's factor solved.
And there are some very important thingsthat actually to know about it.
And some, if you hadn't seen these
stuff before, you'll learn somethingshortly and it'll
be quite real I mean, it'll be somethingyou should know, and it's not expected.
So the idea is this, there's two
steps, there's factor and there's solveAnd factor
says you take your coefficent matrix andyou
factor it as a product of simple matrices.
By the way, you've seen tons of thingslike that.
Like SVD, eigon decomposition, eigenvector, you
know, you've seen lots of these things.
Right?
QR factorization.
Alright, so what you do is you factor itas a product of simple matrices.
Sometimes it's two, it just it could be asfew
as two and it could be as many as four orfive, right.
Actually sometimes it can be a variablenumber where some of those things
have exotic data structures, like it couldbe these products of reflections, right.
OK.
So well, one example is you can actuallyunderstand the fft this way.
I don't knwo if people know.
If you know what the fft is, great andI'll say what it is.
it's a specific matrix.It's the dfd
matrix, the discreet
[UNKNOWN]
transform.
It's a Vandermon matrix with a bunch ofcomplex numbeers in it.
Do you know what I'm talking about?
Fine.It basically calculates a discreet
[UNKNOWN]
transform.
And there's an interesting thing where youtake an order n, matrix like that, and you
can factor it into a product of log
and matrices, which have very specificsparsity patterns.
Okay, and it turns out, if you do this andthen
exploit some of the structure, you get anN log N factorization.
Okay, so, alright, let, let's move on.
So, these are things like, you know, diagonal, upper, lower, triangular,
and so on.And then I mean this is totally obvious.
>> What is A inverse B?
Well, I mean, the inverse of a product isjust the product of the inverse.
It's in the other order, like this.And then you parse it.
This is extremely important this way,remembering
that you're never computing any of theseinverse.
Right?You're, you parse it this way.
Right?
So, you should actually think of that as amethod.
A1 inverse that you call on B.so.
That, that's what happened.
So, it, you just un, you just, you unwind,right?
So you write A as a product, calc,calculate
inverse B, you call the inverse methodsfrom backwards, right?
And you just unwind everything. So it looks like, that's the idea, right?
So, so that's the idea. now this immediately has some.
Pretty cool implications. It says the following,
and this, this is like very important.
It says that suppose you need to sol, hadsolve a bunch
of equations with the same coefficient
matrix but with different right handsides.
Okay, that comes up in tons of places,right?
I mean Void, tons.Right?
so suppose that's your problem.Right?
And what you, by the way, you'll see tonsof cases where this comes up.
Okay.Then, what it says is, you factor
the matrix.
That typically costs more we'll talk aboutthat
in a minute, then the, then the solve.
Right?
so you factor it and you cash thefactorization.
Now what that means is it says that if youwant, if you have, if you want to
do multiple right hand sides you get onefactorization, here, and m solves.
Okay?
So, for the most common methods by theway,
solving linear equations, the factor plusn cubed, the factor steps.
And the back solves cost n squared thesolve steps.
Okay, and now I'm going to ask you aquestion.
To solve 1000 linear equations with 1000variables on
that machine is, I'll just make up anumber.
We'll fix it as, I like your number.30 miliseconds.
It's 30 milliseconds,
and I want to know how much time does ittake to solve two equations.
AX equals B1, AX equals B2, with the sameA.
It's 30 milliseconds to do one. How many is it to do two?
By the way, the naive method just calls,solve AB twice.
And the answer is, what did we say, 30millisecond?
60 milliseconds.
And what is the correct answer?30.
It's studen, is 30.
OK?All right?
How many
does it.So, by the way, that's pretty weird.
Everyone got that, that the cost ofsolving two linear equations, if you,
if you, if you cash the factorization, isactually the same as one.
And I don't think this is obvious,frankly.
it is not obvious unless you knew this,right?
It, it's an absolutely real phenomenon andit's, it's the key to tons of things.
Okay?Because of things being sort of either
viable or not.Right?
how much would it cost me to solve, let'ssay, 25 linear equations.
25 different righthand side same A on thatlaptop.
I'm doing 25 times more work. And time is the same.
The answer is yes because thefactorization dominates everything.
Once you factor you can it simply.
If whatever your background is, you didn'tknow it, this
can actually have a huge impact on allsorts of stuff
So let's look at some factorizations.
By the way, I'm not going to
cover factorizations here so that's anentirely different
topic and they're not that hard actually
to, to work out how these factorizationsgo.
It's a whole world in itself.
But this just so you know at the highestlevel.
We can't go into the details offactorizations.
So.Here's the most famous one.
Is LU factorization and by the
way it's completely equivalent to gaussianelimination.
It says any nonsingular matrix you canfactor
as PLU.
Where P is a permutation matrix, L islower triangular,
U is upper triangular and obviously L andU have
to have non-zeros on their diagonals causeif there was
zero on the diagonal of L or U, what wouldhappen?
Yeah.It's non-invertible.
The determinant is zero and so you have noright solving Ax equals b, okay?
So, okay.
So, and this factorization,
we're not going to go into it although itis
just Gaussian elimination, the cost isorder n cubed, FLOPs.
Right?
So, that, that's how that works.
Okay?And so, here it is.
Here's how you solve a set of linearequations by LU factorizations.
So, you're given A and B.
The first thing you do is you factor A asPLU.
That's N cubed FLOPs.
you then, then how do you solve? Well, you go backwards, right?
So,
you solve PZ one equals B, that's like
zero flops or whatever, I mean, whateveryou like.
It's, you premuting the entries of B.
Right, then you do forward substitutionback substitution, and
you can see this is like two N squared.
Right.
And so, I mean the numbers don't reallymatter
here but the point is n is typically bigenough.
N would have to be tiny, before like 10 orsomething, then all of this is off anyway.
But n would have
to be, so if n is any reasonable number10, 100, 1000
or something like that, then the n cubedthing is bigger, right?
And so that's the idea.
And so it's basically about n cubed,right?
It's, it's, it's the, it's, it's thecomplexity, okay?
So that, that, and that's LUfactorization.
the back solve, notice it's order nsquared.
So, in fact, you would say the following.the, let's see,
how shall I put it.
Okay, here's a good way to conceptualizeit.
one way to think of solving a x equals
b for different right hand sides, issomething like this.
The first time you solve it, it costs youn cubed.
Thereafter At, after you ca, you factor Aand cash it,
you get a discount by an or, by a factorn.
Everybody got that?
Right?
So, so the first one might take you half asecond and then the next
one is going to take you submillisecond,right?
I mean, so these are actually veryimportant things to know about.
Okay.
now, sparse LU factorization. this is quite interesting.
And I'm not going to get into the details,
but this is actually very good stuff toknow about.
if any, if you have a sparse matrix A, themost common
factorization method is to factor it intoa product of four matrices.
and you know by the way they have othernames for these you know p1 and p
2 were sometimes called you know row andcolumn
permutations or something but that is whatthey do.
That is you factor it into a product offour things a permination LUP okay.
Now the fact is mathematically you
only need one of those permutations toensure the existence of a factorization.
Right.
That's the mathematically because itfollows from the previous one.
Right, the reason you'd use two and reasonyou would use a permutation both
on the left and on the right is that the Lthat's actually happening here is this.
As P1 transpose A P2 transpose
Right is L U.
So what you're really doing is permutingthe rows and columns
of A, and then calling an unpermuted L Ufactorization on it.
Okay?
Now, by the way, if you choose P1
and P2 wrong, L U factorization doesn'texist, right?
So let's assume you've chosen themcorrectly and they work.
It turns out that A is sparse.
The sparsity pattern of l and u
is, depends on the permutations.
Okay so if you pick a good permutation land u are quite sparse.
And that, of course, is going to affectyour solve time.
Because, your what, your solve time. We've already looked at that.
Back and forward substitutions.
Basically n and z.On L, and then Z on U.
Alright.
That's the number of flops.
So the point is, if I pick 1 permutationP1 and P2, so that, L has
a milion non zeroes and I pick another one
so that L has a hundred million nonzeroes.
Right.
It's not just a memory issue but also arun time.
Everybody see this?
So okay, this is actually quitefascinating It turns out
that when you do an LU factorization on amatrix, right?
Then the permutation can have a hugeeffect on the sparsity of the result.
How do you choose P1 and P2?
Well, you can ask lots of questions.
For example, obvious question would be,how
do you choose row and column permutationsof
a sparse maker so that it's, it's
LU factorization has the fewest numbers ofnon-zeros.
It's a totally natural question and it'smp heart.
Okay?So, fine.
so what do people do?
Well, it's, it means everything, these areall heuristics.
Okay?
And, actually, they're more thanheuristics, right?
it's, they're even,
they're very stringent requirements,right?
Because what people do who worked onsparse matrices is the following.
If you want to, so what happens is, something like this.
Well, you choose these entries usingvarious methods, right?
One of the methods depends on the sparsitypattern,
it also depends on the values of A, right?
As you're going along you might choose,the P
is actually built up while you're doingthe factorization.
Okay?
So that, that that's how this works.
And, and I'm not going to go into anydetails, but that's sort of the idea.
And it's, it's worth knowing about this.
Okay.And, what's the cost?
Well it's usually, way less than n cubed.
Right.
By the way, if it gets even close to ncubed It means you have
implemented an extremely inefficientmethod when you should
have just called a, a dense LU packroutine,
right, which would've blocked everythingnicely for you.
That's, that's all it means, right?
But, but this is a big deal. so you know, if you're lucky this is
>> This is actually I'd, a, a, a verybig
give, by the way, everyone here is abeneficiary of this.
I can tell you that right now.
so, and the reason is very simple.CVX, for example, calls SDPT3 or SEDUMI.
All of these things use sparse matrixmethods to actually do the heavy lifting.
That's if you profile any of thesesolvers, that's all they're doing.
The only thing they're doing is solvinglinear equations.
So, there's nothing else.
Linear, And they're literally solving a xequals b.
Now a is symmetric, but that's anotherstory, it doesn't matter.
Same idea goes, right.
And they're all spa, and, and they'reexploiting sparsity very heavily.
Okay, now they, what's at, there's alittle, I mean it's a bit funny, right,
because.
You know, if someone comes along over yourshoulder, when you're
plot feeling with CVX, and somebody says,what are you doing.
Well, I'm doing convex optimization, andthey go, well, you know, how
do you know it works, and you go, that'sthe whole point.
It's just like, it's, it's tractable.
I mean, we can do that, I mean, you know,
had a whole course, and that, that was thewhole point.
There's other problems you can't solve, but these you can really solve.
The, as a practicle matter the only reasonyou can solve them, as a practicle matter,
is because the euristics of choosing thepermentations
for solving linear systems of equationsare so good.
Everybody got that?
Don't, but do not tell anyone outside thisclass that.
Right?This is a secret we keep together here.
Right?
So, if anyone asks you, what, why doesthat work?
Is that heuristic?
You say absolutely not.
This is not a heuristic. and in theory it's not,
right.
If, the worse thing that can happen isyour
permutation fails spectacularly and youavert to n cubed.
That is not good.
[INAUDIBLE]
10,000 or 50,000.So, but now you know.
This is why, roughly speaking, you cansolve equations.
10,000 variables, 100,000 variables, amillion variables.
Million by million.You can, kind of, solve, sparse equations.
And, so, the correct way to say it is ifthe gods, god or gods, who control The
[UNKNOWN]
of permutations used in sparsefactorization are smiling on you, and
that depends on the offerings you havemade to them likely.
If they are smiling on you then you cansolve equations that would just otherwise
be completely and totally out of question,everybody
got this right but there are some sparsity
[UNKNOWN]
Where you don't have to rely on
the gods of heuristic, sparse permutationselection algorithms.
And we'll talk about some of, well I'lltell you one right now.
If the matrix is banded, you don't need tomake any offerings at all.
That's, it's really going to work.
Always, m'kay.
But you know for general sparse matricesthis is not true, okay, so.
Now this is very
good stuff to know about even if you, ifyou're more interested
in sort of more theoretical stuff, this isvery good stuff to know about.
So now, well we're going to encountermore, I mean what we're going to encounter
are some things that are generally not
just square, but you actually have sinestructure.
There's the Cholesky factorization.
If you have a positive definite matrix, A,then you can factor it as it LL
transfer \g, but it's not.
What it says is, it's an LU factorization,with L = U.
Well kind of duh, or aesthetics tells youthat right.
because what's a matrix that's symmetric.
Well it's transposes the same, or you swaprows and columns, the same thing.
So, of course, the U and the L, sh, that'sall this is, okay.
And L is lower triangular And in
fact, generally people require just tonormalize that.
That the entries are positive.
One very interesting fact about Choleskyfactorization, is a numerical one.
It actually is stable.
Again we're not really talking about that.
It means you do not have to permute, you,you
can say immediately, that if you calculatethat in double precision.
You aren't, you're not going to have any,any trouble, right.
That your er, you can bind the error in L.
And stuff like that.
As a, as a practical matter, it justworks.
In this case, and
the cost is 1/3 n cubed. And, you know?
Again, the difference between 2/3 and 1/3
is, like, zippo in, in current things,right?
Where cache, cache misses is way moreimportant than some factor of two.
But, you know?
You get the, rough.
[INAUDIBLE]
, I mean, this, you save because.
Well, you don't have to calculateseparately in an L and a U, right?
You're just calculating one L.
Okay.
So, how do you solve things by Choleskyfactorization.
by the way, this would come up in leastsquares, right?
That's your most, that's your most
basic convex optimization problem, leastsquares.
You can do it by Cholesky factorization.
So.
Here, you'd factor A, as l l transpose,
and you do forward and backwardsubstitution, and
the cost is about 1 3rd n, cubed.
Okay, so that's, that's the, that's theidea, that's Cholesky factorization.
Now sparse Cholesky factorization.
it looks just like sparce LU, but there'sactually a
very critical difference, and I'm going tosay what it is.
And you probably don't need to know it at
this level, but these are good things toknow about.
So here, here it is.
as before, there's P one and P 2, you dopre- and post-multiply.
Or, if you like, you can say that whatyou, another way to
write this Is, is to say just to think ofit this way
first you permute A then call a Choleskyfactorization on their perumuted version
right so that, that, that's a a good wayto think of it.
Oh and by the way.
Everybody has a different thing, you know.
Sometimes it's lower triangular, sometimesit's upper triangular.
The permutation either goes like
here, or there.
So you always, anytime you actually look, use something that does a factorization.
You have to look at what, exactly whatthey mean by it.
And it'll probably be different from thesenotes.
Right, so, but.
That's, this is too bad. There's no standards.
Okay, now here you always have a Choleskyfactorization.
But the, the permutation can have a hugeeffect on the sparsity of L, and therefore
on the memory required to carry this out, and second, and of course, it
also has, since the number of nonzeros inL is basically the flop count.
That's going to affect the run time.
Okay.
Now, unlike the sparse LU factorizationwhere you're actually choosing
P1 and P2 as you go, as you factor thematrix.
Right?
Because you actually looking at entriesand
you want to, you don't, want to avoidnumbers.
They're too small.
If you wanted to, if you're going todivide by something.
Here it turns out that P can actually bechosen
entirely on the basis of the sparsitypattern of A.
Okay?
So, in other words, all that mat-, youtake the sparsity
pattern of A, you don't even have to knowwhat the
entries, the numbers are going to be, Right, if you can make
a graph out of it that's a very commonthing, and
then you chose an elimination orderingright but you do this once right.
So, that says if you are going to dosomething like
solve ax equals b with different values ofA multiple times.
When A is sparse definite, you canactually cash the permutation.
I don't know if you got that but this willbe,
this will, the usefulness will be clear inabout two weeks.
You can cash the factorization.
and that's the idea.
anyway, so the cost is usually you know,way
less than and and all that kind of stuff.
So that's, that's how that works.okay,
another factorization we going to see isthe so called LDL transpose factorization
and that's this if you have a symmetricmatrix not necessarily not necessarily a
positive definite you can factor it as LDLtranspose right and D.
P is permutation matrix.
l is lower triangular unit, so you have li
equal 0, d is block diagonal, and it'seither got
1 by 1 blocks or 2 by 2 to blocks, right?
So and you get sort of the same thing,right?
So this would be kind of one, this would
be one of the methods of, of sort ofchoice.
Oh, and I think I, I can give you a roughidea of what
>> What backslash doesn't MatLab.
I mean, just roughly. It's something like this.
When you type A backslash B, I think itlooks at the diagonal of, of A.
And if they're positive,
or it looks, it quickly looks.
It compares some entries quickly to seeif, if
it's possible that A is positive symmetricpositive definite.
Right?
So, it guesses that.Then, it attempts the Cholesky factor.
If it finds that it is, if it, if there'ssome evidence that A is symmetric,
and positive definite, like, for example,it has
a negative on the diagonal, it doesn't dothis.
It's all over.
It's not positive definite.
So, it attempts to do a Choloesky
factorization because it's numericallystable and fast.
Okay?
But, if it turns out A wasn't, then theCholesky fails.
Then, it falls back and does somethinglike an LU or something like that.
I mean this is roughly what it does.
Right?So, okay.
And it won't a sparse factorization unless
A is like explicitly stored flagged assparse.
Okay.Now we're going to talk about something.
It's about an elimination method, but itactually ties to a whole bunch of stuff.
And it turns out it's actually kind of the
same as selecting the permutation, butwe'll get to that.
so what we'll do is let's just take xequals b and block it.
So, all we're going to do is
we're going to, we're going to write it asblocks.
we're going to partition x that'llpartition b, and we're going to write it
that way.So, it's just two sets of equations.
It's a1 1x 1 plus a1 2x 2 equals b1, and
then a2 1x 1 plus a2 2x 2 equals b2,right?
So that's, that's it. On this block, you know.
And let's see what happens.
What you do is from the first equation,basically, is a one
one, x one plus a one two, x two equals bone.
And what you do is assumingn a one one isnon singular, right?
You put this
over there then premultiply the a one oneand
that, that says that you can write x one.
As a function of X2 and you get this,okay.
And so the idea is now we're going
to eliminate X1 from the system ofequations, right.
So, all we do, is we plug this value of X1into this, right.
by the way, this is, this is actuallyjust, sort of, like you know, this,
this substitution stuff.
And with, but with matrices, blocks,right?
So, what happens is, if you plug thisequation into that, right.
You get A 2 1 times this thing, plus A 2 2times x 2 equals b 2.
And when the smoke clears, that's a set ofequations for x 2.
And you get.
This equation right here okay.
So this just you've just elim and you
would say things like you have eliminatedX1 okay.
This matrix you have seen before that'sthe shirt compliment
okay and you've seen it in any legitimateclass you've taken.
Because a legitimate class uses linearalgebra, and then,
you can't use linear algebra and not seethat, right?
So you take probability class, that'sconditional, that's conditional, you
know, covariance or some, that's
the conditional variance, covariancematrix, right?
You take other things.
Take circuits, mechanics, anything,physics.
This comes up.
Okay, so that's the sure compliment.
so what it says is that when you eleminatevariable, a block
or variables, or a block of variables froma set of linear equations.
You're left with a set, with a
linear equation that involves the surecompliment, okay?
So actually, that's, I mean that's goodthing to know about, okay.
So, by the way you would form this matrix, you'd solve X2 and then
you would recover X1 using that formula,okay.
And so, here's a completely generic,completely
generic algorithm for solving this blockequation.
And it says this It says, yeah I
mean you're going to form, you're going toform this matrix.
that's the sure his thing compliment, it'scalled S here.
Right, so you, you form this thing, a 1 1inverse a 1 2.
Now, by the way, that can be formed lotsof different Wait oh,
let me ask that.how would you if there were no structure
in a11 and a11 was let's say a 1000 by1000, and suppose a12 was 1000 by 10.
Please tell me how you would form a 11inverse, a 12.
What would you do?
So A 1 1.Thousand by thousand.
Inverse times,
A 1 2 is a thousand by 10.
How do you, how do you calculate thatproduct?
It li-, it's a product.
How do you do it? A 1 1 cash.
And do 10 back solves.Done.
Perfect.Okay?
So.All right.
So, but, you just write, people just writestuff like form A11 inverse A12.
Okay.
And assuming you know what you're doing.Okay.
so, you form this, your compliment.
then you solve this, your complimentequation, and
then, you use this equation to reconstructX1.
Okay?
So, that's the method.All right?
And, you know let's, let's, let's figureout how much it costs.
And the answer is, it's F plus N2S.
That's the cost of factoring A11.
And S is the solve step.Right?
these are the back solves.Right?
As you, as you said, that there's onefactorization and some solves.
step 2 cost this.
And what we're going to do here is we'rejust going to
assume on the second step as dense.
I mean, because you could work it out, when it's not Right and then the
the the last step is a something that'sagain we're going to just assume its
dense generic nothing and we'll just do
we'll just solve that equation and that'sgoing to
be a N, N, N two cubed right so the totalthing looks something like this.
Now.
Here's what happens.If, if you apply this method now, to a
matrix that has no structure whatsoever, then f equals this, and s equals that.
And you plug it into this formula, and youget back, actually, exactly.
The fact you get, there's no savings whatso ever okay?
So that was a completely silly thing todo.
I mean, it doesn't matter that's how youdo it right.
Actually, it's interesting because itshows you,
it shows you what blocking is right?
because there is, this is blocking 2 by 2now you can imagine what blocking
10 by 10 looks like, it's just this. Okay, so that's it.
Here is where it's useful.
>> When A11, has structure, right?So, the factor, with the factor cost and
the solve cost on A11 is low, then, thisdelivers, this delivers big savings.
Here's an example. Suppose that A11 was diagonal.
Right?That means the factor cost
is zero and the solve cost is n1. That's just, just divide.
Right?
Then the number of flops is this, and thisis absolutely amazing,
so this is that if you see a set of linearequations
that look like this.Right?
So, It looks like that.
I'm going to, I'm going to draw the,non-zeroes.
Right?
Oh sometimes people call this a spydiagram,
I have absolutely no idea why, but anyway,whatever.
I, it doesn't matter.
I'm just showing you the non-zeroes.
Okay?
Now there's many other ways to understandthis,
and and I'll, I'll show A, a couple and
the once the way they are purely
mathematical you probably have seenbefore, OK.
So and also introduces a very interestingidea so you what that is.
So lets , lets takes this one suppose youwanted to solve A plus BC and
the idea here is that B is like you knowit, its so this is low rank.
Looks like that.
And here's, and here's a. and let's assume a has structure.
So let's make a banded.You know?
That's, 'cuz it's pretty, and it hasstructure.
So, there's your, your sparsity pattern.
That's a plus bc, okay?
it could be anything else. It could be sparse.
It could be.
Forty eight, you know, DFT, it could beanything, you know,
it could be something that you have a fasttransformer for.
It could be some wave lit thing orsomething, I don't care, right?
And then so the people would call this islike structured plus low rank.
That's what this is.
So it could be sparse plus low rank matrixor something.
Okay.
Now the dif, the sparsity pattern of Aplus BC is generally full.
Severy, you know, cause, you take anouter,
you just take B and C to be columns.
Or, you know, B is a column and C is a rowvector.
Then it's just an outer product, it spraysentrys all over the place.
Right, so, so then, a generic, a genericmethod there will be no sparsity, nothing.
So, how you can solve these is veryclever, and
in fact the first step, there's a name forit.
It's actually called unelimination.
So elimination means you take yourvariables
and you look at one block and you
say, I'm going to express that in terms ofthe others and then get rid of it.
Unelimination means that you're going tointroduce a new, a new one.
So what you do is you look at (A + BC)x =
b, you write that as A + By.>> Right?
Equals B.Y is actually CX.
Everybody got that?
That's just, of course, they're the samething.
They're identical.
But, you rewrite these equations likethat.
Right?
So, there's, there's, and , and, here, so, you've, i mean you've done unelimiation.
You've introduced new variables and newequations.
Everybody got this?
And you see this.
Now, you're not trained yet but you willbe.
And it, once you're trained, what youreye,
what you're searching for is the
following: easily inverted matrices on thediagonal.
And your eye should, once you're trained, go immediately to that minus I.
You will go, it will go immediately there.
Because in the world of easily invertedmatrices,
I is up there and -I is up there.
Okay?Everybody cool?
Alright,
now.
I need, see, I just have this fantastic, Ijust learned about block elimination.
There's a, I'm going to elimi- I eliminatethat block.
Now if you eliminate that block you'reeliminating y and where do you end up?
Right back here.
Everybody got it?
So, well duh, I mean you un eliminated and
then eliminated the thing you had justintroduced right?
OK.
Fine.So the key, is to, is to control.
>> Your urge.
Primal matrix structure factorizationurge, to eliminate that minus I.
And instead, you take the system and youeliminate A.
And, in fact, that's not a bad ideabecause
the assumption here is that A is easy toinvert.
Right?
So, A is 48 transform or some wavelifting.
I don't know what it is. It doesn't matter what it is.
Right?
So, that's what it is.OK.
If you eliminate A, so, if youuneliminate.
If you go from this.To this, and then you eliminate a.
Not the minus i.
You get this equation.
And you will recognize that as the Schurcomplement here, right?
It's like, you now, this thing minus thattimes the inverse of that times that.
OK.Or other way around.
Some, that's that's Schur complement.
OK?And you get this.
Now, then you solve for ax, and you can,once you get y,
you get x, okay?
And you can actually write this out, just
purely mathematically, and it comes out tobe this.
A plus b, c inverse, is a inverse minus, and then this big old thing over here.
Okay, it's a very thing, it's got allsorts of names.
matrix inversion lemma I think is a niceone, and you'll just
do a quick example and then, and thenwe'll quit for today.
so if A were diagonal and B and C weredebts then
that's a beautiful example, right?
So that's called a diagonal plus low rankMatrix, right?
Perfect, perfect example of thisstructure, right?
then it says you can actually solve it inhere
if you, you do the dumb method it's like ncubed.
If you'd used this method which is anun-elimination and an elimination it
actually ends up it goes from cubic in nto linear in n.
OK?
So my claim is, you just, these are not
small differences, right, these are nottheoretical things, right?
The difference between knowing this andnot knowing it means you can look
at something, and it could be a million bymillion set of equations.
But if it's diagonal plus low right, nowyou know.
I can solve that in probably in likesub-millisecond.
But you have to know what you're doing.