Column Vectors vs. Row Vectors(steve.hollasch.net)
steve.hollasch.net
Column Vectors vs. Row Vectors
http://steve.hollasch.net/cgindex/math/matrix/column-vec.html
78 comments
I was with you until you said bra and ket didn't do it for you. Its right there in the name bra(c)ket. <bra|, |ket>. Logically equiv to any other notation, but impossible to fudge up. When they are multiplied in a way that lines up like <bra|ket>, you got an inner product (scalar measure of overlap). When they are multiplied wrong way around, you got a ketbra |i><j| (which is basically a neat way to represent a matrix element at indices i and j). Is the bra the row and the ket the column or vice versa? Don't matter. Just keep things facing whichever way they were originally facing and nothing will go wrong. In fact, you don't even need to think in terms of columns at all if you don't want. That's it. That's the whole system.
For the engineers who never seen this, |i> |j> |k> work exactly like unit vectors i,j,k, but you aren't limited to the usual 3 dimensions and 3 letters. Writing |x> in a ket like this just means "I am using x as a label for a degree of freedom in some vector space, not as the literal value x". This lets you do cool shit, like if you got n dimensions you can write a vector as sum v_i |i> ranging from 0 to n. You can also do super useful stuff like merging two vector spaces together to represent a larger system (aka a tensor product). Eg if I have a vector v = sum v_i |i> representing one subsystem and another vector w = sum w_j |j> representing some other subsystem, I can construct the vector u = sum v_i w_j |i, j> representing the whole system (sum is over both i and j).
(I know you already learned all this, I wrote it for everyone else's benefit).
Unrelated, given the webs origin at cern, HTML should have used bras as start tags and kets as end tags. <HTML| <head|head> <body| ... |body> |HTML>.
For the engineers who never seen this, |i> |j> |k> work exactly like unit vectors i,j,k, but you aren't limited to the usual 3 dimensions and 3 letters. Writing |x> in a ket like this just means "I am using x as a label for a degree of freedom in some vector space, not as the literal value x". This lets you do cool shit, like if you got n dimensions you can write a vector as sum v_i |i> ranging from 0 to n. You can also do super useful stuff like merging two vector spaces together to represent a larger system (aka a tensor product). Eg if I have a vector v = sum v_i |i> representing one subsystem and another vector w = sum w_j |j> representing some other subsystem, I can construct the vector u = sum v_i w_j |i, j> representing the whole system (sum is over both i and j).
(I know you already learned all this, I wrote it for everyone else's benefit).
Unrelated, given the webs origin at cern, HTML should have used bras as start tags and kets as end tags. <HTML| <head|head> <body| ... |body> |HTML>.
Fun anecdote: when I was a teenager learning LaTeX, I'd got in the habit of google image searching the symbols – "Latex infinity" would turn up a pic of an infinity symbol with "\infty" below it; "latex contour integral symbol" would show "\oint", and so on.
As a nineteen year old first-year university student, I genuinely typed "latex bra" into google image search, pressed enter, and then had a second or two of shock as my brain was genuinely surprised to not see quantum mechanics in front of me with a nifty little "\langle" involved, but rather lots of women in decidedly unique looking underwear...
As a nineteen year old first-year university student, I genuinely typed "latex bra" into google image search, pressed enter, and then had a second or two of shock as my brain was genuinely surprised to not see quantum mechanics in front of me with a nifty little "\langle" involved, but rather lots of women in decidedly unique looking underwear...
I learned that "garbage collection" is not just a memory management technique but also refers to big trucks removing the content of trash cans when I took a "nerd test".
One of the architects of the Pentium Pro, Intel’s first out-of-order execution engine, told of walking up to a Coke machine with an “out of order” sign on it — “Oh, out-of-order has other meanings.”
The vending machine at intel drops a can of coke. A CPU architect walking by picks it up then inserts his money. The machine is out of order.
I am actually a fan of bra-ket notation. But in my mind, <i|j> is the amplitude for going from state i -> j. My point is, if you just have one state that you want to talk about, which one do you write down <a| or |a>? It doesn't matter. You only write the brakets when you want to calculate transition amplitudes/inner products, or cross products.
I agree about just keeping them the way they are originally facing and nothing will go wrong. The nice thing is that you can mechanically rearrange things almost like regular algebra, i.e. draw out a common factor, and don't have to think what these objects are (a vector? a dual vector? a functional?).
Same with four-vector notation / Einstein summation. As long as you are following the rules for raising and lowering indices, you cannot form "wrong" products or mess up the transformation properties.
Said in another way: as a (former) experimentalist, I don't care if the vector is p^mu or p_mu. All I care about is the product p^mu p_mu.
I agree about just keeping them the way they are originally facing and nothing will go wrong. The nice thing is that you can mechanically rearrange things almost like regular algebra, i.e. draw out a common factor, and don't have to think what these objects are (a vector? a dual vector? a functional?).
Same with four-vector notation / Einstein summation. As long as you are following the rules for raising and lowering indices, you cannot form "wrong" products or mess up the transformation properties.
Said in another way: as a (former) experimentalist, I don't care if the vector is p^mu or p_mu. All I care about is the product p^mu p_mu.
If you define u to be a column vector, and u' to be its (conjugate) transpose, and define matrix multiplication as usual, then u'v is a dot product and uv' is an outer product. What you mentioned about the standard basis vectors in R^3 works for standard basis vectors e_i in R^n using this notation; e.g., e_i e_j' is an all zero matrix with a one in the ith row and jth column.
If you have a hard time remembering whether u should be a column vector, you can simply identify things that look like u with vectors (|ket>), and things that look like u' with linear functionals (<bra|). Then u'u is just a dot product.
There are other notations, too... Sometimes (u, v) and u \cdot v are inner products, sometimes u \otimes v is an outer product.
Generally linear algebra-style notation is terser and has less line noise. IMO, some manipulations are simpler to think about using this notation than with the bra-ket notation, especially when you're doing numerical optimization or numerical linear algebra.
Maybe another reason to prefer it is that modern notation for operator algebra follows this notation and not the bra-ket notation (although I think bra-ket notation was invented for this purpose...). If L is a linear operator, L^* is its adjoint, Lu is the application of L to u, etc...
Of course, what notation you should use is determined by what community you belong to. If you're publishing papers on numerical linear algebra, it would be very strange to decide to use the bra-ket notation instead. The purpose of notation is communication. The goal is to communicate ideas to other researchers. It's good to develop the ability to be flexible with notation and go along with what other people are doing.
If you have a hard time remembering whether u should be a column vector, you can simply identify things that look like u with vectors (|ket>), and things that look like u' with linear functionals (<bra|). Then u'u is just a dot product.
There are other notations, too... Sometimes (u, v) and u \cdot v are inner products, sometimes u \otimes v is an outer product.
Generally linear algebra-style notation is terser and has less line noise. IMO, some manipulations are simpler to think about using this notation than with the bra-ket notation, especially when you're doing numerical optimization or numerical linear algebra.
Maybe another reason to prefer it is that modern notation for operator algebra follows this notation and not the bra-ket notation (although I think bra-ket notation was invented for this purpose...). If L is a linear operator, L^* is its adjoint, Lu is the application of L to u, etc...
Of course, what notation you should use is determined by what community you belong to. If you're publishing papers on numerical linear algebra, it would be very strange to decide to use the bra-ket notation instead. The purpose of notation is communication. The goal is to communicate ideas to other researchers. It's good to develop the ability to be flexible with notation and go along with what other people are doing.
Unfortunately this physics-guy/spherical-cow thinking is really common in people doing linear algebra on computers and the results are poor performance. Abstractions are great for general programming, but when you are doing tight-loop math it's really not applicable. I've seen this a lot recently with "clever" NDArray style math libraries (not to mention the utility of more than 2 dimensions is .. very limited)
In my not-super-extensive experience, if you're programming and dealing with linear algebra problems - and you're past the MATLAB prototyping stage - then I really suggest using the BLAS/LAPACK API. You're not gunna beat decades of nerds programming ballistic missiles on punchcards. They're kinda weird and unwieldy and don't match to textbook math - but they've had a ton of thought put into them. And are made with as few footguns as possible. At least re-implementing an algo with them a few times is really educational. You'll see that it forces you to think about your memory layout and you'll realize that the back of napkin complexity calculations are actually tricky to massage into a good algorithms that use memory well. The final result with something like the Intel MKL will be way way faster than anything your abstraction can achieve
(Unfortunately I haven't found a good primer on using the BLAS well)
In my not-super-extensive experience, if you're programming and dealing with linear algebra problems - and you're past the MATLAB prototyping stage - then I really suggest using the BLAS/LAPACK API. You're not gunna beat decades of nerds programming ballistic missiles on punchcards. They're kinda weird and unwieldy and don't match to textbook math - but they've had a ton of thought put into them. And are made with as few footguns as possible. At least re-implementing an algo with them a few times is really educational. You'll see that it forces you to think about your memory layout and you'll realize that the back of napkin complexity calculations are actually tricky to massage into a good algorithms that use memory well. The final result with something like the Intel MKL will be way way faster than anything your abstraction can achieve
(Unfortunately I haven't found a good primer on using the BLAS well)
MATLAB uses BLAS and LAPACK under the hood. You can write reasonably fast code in MATLAB if most of your time is spent calling out to faster low level implementations.
There are also many different implementations of BLAS, with varying performance. The API itself is not the reason for its performance. Actually, I would argue that it’s a pretty crufty API which should probably be done away with at this point and replaced with something less awkward. I don’t know how many times I’ve had to remind myself what a “leading dimension” was…
A good example is BLIS, which does provide the BLAS API but has a more modern “object-style” C API which is significantly easier to work with than the BLAS API with no performance hit.
Also, depending on what you’re doing, abstractions can be very helpful indeed to keep around. Even at the BLAS level, groups of bits are thought of as singles, doubles, real or complex… a matrix has a size and shape… etc. Having this information is useful. Linear algebra is full of type information which can be used to dispatch different algorithms.
Lots of other alternatives to BLAS which are a bit higher level but still very useful: Eigen and Armadillo in C++, Julia, etc.
Anyway, I would say (based on my actual, significant experience ;-) ) people are beating the guys with punch cards all the time, no reason to stay in the 70s…
There are also many different implementations of BLAS, with varying performance. The API itself is not the reason for its performance. Actually, I would argue that it’s a pretty crufty API which should probably be done away with at this point and replaced with something less awkward. I don’t know how many times I’ve had to remind myself what a “leading dimension” was…
A good example is BLIS, which does provide the BLAS API but has a more modern “object-style” C API which is significantly easier to work with than the BLAS API with no performance hit.
Also, depending on what you’re doing, abstractions can be very helpful indeed to keep around. Even at the BLAS level, groups of bits are thought of as singles, doubles, real or complex… a matrix has a size and shape… etc. Having this information is useful. Linear algebra is full of type information which can be used to dispatch different algorithms.
Lots of other alternatives to BLAS which are a bit higher level but still very useful: Eigen and Armadillo in C++, Julia, etc.
Anyway, I would say (based on my actual, significant experience ;-) ) people are beating the guys with punch cards all the time, no reason to stay in the 70s…
> MATLAB uses BLAS and LAPACK under the hood. You can write reasonably fast code in MATLAB if most of your time is spent calling out to faster low level implementations.
Many years ago I was implementing a statistical model in R that involved a lot of matrix multiplication (the model included a non-negative matrix factorization step). Fed up with how slow it was in R, I quickly reimplemented it in MATLAB--literally copy-pasted the code and changed things like operators and function names--for an immediate 5x speedup.
Many years ago I was implementing a statistical model in R that involved a lot of matrix multiplication (the model included a non-negative matrix factorization step). Fed up with how slow it was in R, I quickly reimplemented it in MATLAB--literally copy-pasted the code and changed things like operators and function names--for an immediate 5x speedup.
"The API itself is not the reason for its performance"
It eliminates many performance mistakes by making it so that when you do something dumb it's generally immediately apparent. This mostly boils down to nonsensical/slow operations being usually impossible and any memory copying being explicit.
I can only tell you my own experience. Going through the steps of transforming a few "textbook algorithms" into BLAS made me think about the problems much more clearly and in a way I would have missed using a higher abstraction.
All the examples you give (Matlab, eigen, armadillo, Julia) allow you to really easily write really bad code (blas under the hood is kinda irrelevant). But you're totally right that they're useful. If you want to call out to code that does an SVD or something simple that exists in a library then they're generally just fine and you can't mess it up too much.
I use BLAS with a repl from Clojure (thanks to neanderthal) so it's very painless. Sounds like you're using it in a compile/run loop from C/C++ which sounds like very much not fun.
It eliminates many performance mistakes by making it so that when you do something dumb it's generally immediately apparent. This mostly boils down to nonsensical/slow operations being usually impossible and any memory copying being explicit.
I can only tell you my own experience. Going through the steps of transforming a few "textbook algorithms" into BLAS made me think about the problems much more clearly and in a way I would have missed using a higher abstraction.
All the examples you give (Matlab, eigen, armadillo, Julia) allow you to really easily write really bad code (blas under the hood is kinda irrelevant). But you're totally right that they're useful. If you want to call out to code that does an SVD or something simple that exists in a library then they're generally just fine and you can't mess it up too much.
I use BLAS with a repl from Clojure (thanks to neanderthal) so it's very painless. Sounds like you're using it in a compile/run loop from C/C++ which sounds like very much not fun.
Here's BLIS's object API:
https://github.com/flame/blis/blob/master/docs/BLISObjectAPI...
Searching "object" in BLIS's README (https://github.com/flame/blis) to see what they think of it:
"Objects are relatively lightweight structs and passed by address, which helps tame function calling overhead."
"This is API abstracts away properties of vectors and matrices within obj_t structs that can be queried with accessor functions. Many developers and experts prefer this API over the typed API."
In my opinion, this API is a strict improvement over BLAS. I do not think there is any reason to prefer the old BLAS-style API over an improvement like this.
Regarding your own experience, it's great that using BLAS proved to be a valuable learning experience for you. But your argument that the BLAS API is somehow uniquely helpful in terms of learning how to program numerical algorithms efficiently, or that it will help you avoid performance problems, is not true. It is possible to replace the BLAS API with a more modern and intuitive API with the same benefits. To be clear, the benefits here are direct memory management and control of striding and matrix layout, which create opportunities for optimization. There is nothing unique about BLAS in this regard---it's possible to learn these lessons using any of the other listed options if you're paying attention and being systematic.
https://github.com/flame/blis/blob/master/docs/BLISObjectAPI...
Searching "object" in BLIS's README (https://github.com/flame/blis) to see what they think of it:
"Objects are relatively lightweight structs and passed by address, which helps tame function calling overhead."
"This is API abstracts away properties of vectors and matrices within obj_t structs that can be queried with accessor functions. Many developers and experts prefer this API over the typed API."
In my opinion, this API is a strict improvement over BLAS. I do not think there is any reason to prefer the old BLAS-style API over an improvement like this.
Regarding your own experience, it's great that using BLAS proved to be a valuable learning experience for you. But your argument that the BLAS API is somehow uniquely helpful in terms of learning how to program numerical algorithms efficiently, or that it will help you avoid performance problems, is not true. It is possible to replace the BLAS API with a more modern and intuitive API with the same benefits. To be clear, the benefits here are direct memory management and control of striding and matrix layout, which create opportunities for optimization. There is nothing unique about BLAS in this regard---it's possible to learn these lessons using any of the other listed options if you're paying attention and being systematic.
"You're not gunna beat decades of nerds programming ballistic missiles on punchcards" funny and true
On top of that, you almost certainly do not want to instantiate an N x N matrix, because sparsity is an important factor in most real problems.
I imagine that's the right way to do it. Treating a linear map as a table of numbers is too concrete, and there are far too many conventions for how to do it. Abstractly, it's a type-(1,1) tensor, which is naturally isomorphic to a linear map. (I'm slightly abusing how I use the term "naturally isomorphic"). There are many notations for tensors that don't mention indices, including the somewhat surprising Abstract Index Notation by Roger Penrose.
For the sake of completion: A vector is a type-(1,0) tensor. A linear form (or "covector" if you prefer this term) is a type-(0,1) tensor. A bilinear form is a type-(0,2) tensor.
The point of saying the type of a tensor, is to clarify what geometric object an array of numbers represents. A curse (and a blessing) of matrix notation is that it doesn't tell you this. Both a type-(0,2) and type-(1,1) tensor can be represented by matrices, but the first is used for (let's say) dot products, and the second is used for linear transformations.
The only problem with type-(m,n) notation is remembering whether "m" is the contravariant order or the covariant order. It turns out that it's contravariant, which is the opposite to what the naming would suggest. In programming, we would do {cov: m, con: n}, so we wouldn't have to remember.
For the sake of completion: A vector is a type-(1,0) tensor. A linear form (or "covector" if you prefer this term) is a type-(0,1) tensor. A bilinear form is a type-(0,2) tensor.
The point of saying the type of a tensor, is to clarify what geometric object an array of numbers represents. A curse (and a blessing) of matrix notation is that it doesn't tell you this. Both a type-(0,2) and type-(1,1) tensor can be represented by matrices, but the first is used for (let's say) dot products, and the second is used for linear transformations.
The only problem with type-(m,n) notation is remembering whether "m" is the contravariant order or the covariant order. It turns out that it's contravariant, which is the opposite to what the naming would suggest. In programming, we would do {cov: m, con: n}, so we wouldn't have to remember.
If one thinks of vectors as the matrix of some linear map over the scalars (and similarly covectors as linear maps to the scalars), and assuming the usual composition order of maps then it seems there is no choice as to the column/row representation: vectors have to be (n x 1) and covectors have to be (1 x n)
This doesn't say anything regarding their practical storage, but at least conceptually vectors/covectors are a bit more than just a collection of numbers.
This distinction is pretty useful in practice when encoded in types to avoid implicitly mixing things that are conceptually different.
This doesn't say anything regarding their practical storage, but at least conceptually vectors/covectors are a bit more than just a collection of numbers.
This distinction is pretty useful in practice when encoded in types to avoid implicitly mixing things that are conceptually different.
They should not be calles vectors both.
Horizontal vectors are forms (maps from R^n to R). Vertical vectors are vectors (elements of R^n).
Forms applied to vectors are “rows multiplied by columns”, and so: numbers.
Same thing as bras and kets, by the way.
Horizontal vectors are forms (maps from R^n to R). Vertical vectors are vectors (elements of R^n).
Forms applied to vectors are “rows multiplied by columns”, and so: numbers.
Same thing as bras and kets, by the way.
This is a great example of how math is generally weak at types (constantly doing implicit casts and conversions with ambiguous notation, all the way back to Liebniz notation for calculus, and probably further), and the computer scientists add clarity.
It's also an example of how "making it simpler" by dropping types (old Python, Bash, typeless vectors and matrices) ends up making it harder to get work done right when the ideas are non-trivial.
It's also an example of how "making it simpler" by dropping types (old Python, Bash, typeless vectors and matrices) ends up making it harder to get work done right when the ideas are non-trivial.
Well, yes. But it is another way of saying that we think in analogies.
At the same time, forms are also vectors… and vectors are forms of the dual… So it is not that easy.
At the same time, forms are also vectors… and vectors are forms of the dual… So it is not that easy.
If you have two pointy arrows, and want to know if they are perpendicular, then you use the scalar product u⋅v and check if it is 0. There is no difference between u and v in this case. You could say v is a vector, and u a function taking a vector and returning a scalar (a form, a dual vector, ...). But often that is unneccessary artifice. You can just say they are both vectors in the same vector space, which is imbued with an inner product.
Declaring that some objects are forms (duals, covariant vectors, ...) is only neccessary when you want to keep track of transformation properties, for example.
Declaring that some objects are forms (duals, covariant vectors, ...) is only neccessary when you want to keep track of transformation properties, for example.
Abstractly, a vector is neither a row nor a column. An n-dimensional vector is an ordered collection of n coordinates. An m by n matrix is an m-dimensional vector of n-dimensional vectors.
But when you define a matrix in most programming languages you can only really write it like this:
Some languages have fancy syntax that omits the need for braces. But for those that don't, this is clearly an array of 3 arrays of 3 ints. It's only intuitive to think of this as an array of rows given how it's laid out visually.
I also think in that context it would make more sense to index it as mat[row][col], or in other words row major indexing.
Because of that, it also makes sense to store the mat row-major, that is, laid out in memory as [1,2,3,4,5,6,7,8,9] so that mat[row][col] becomes mem[row*3 + col].
This is also how C does it. I think it was a mistake for OpenGL to break with that C convention given how close the relationship between C and OpenGL is. Mathematicians be damned.
But when you define a matrix in most programming languages you can only really write it like this:
mat = [
[1,2,3],
[4,5,6],
[7,8,9]
]
Note that the inner vectors are laid out as rows, and the outer vector (the matrix) is laid out as a column.Some languages have fancy syntax that omits the need for braces. But for those that don't, this is clearly an array of 3 arrays of 3 ints. It's only intuitive to think of this as an array of rows given how it's laid out visually.
I also think in that context it would make more sense to index it as mat[row][col], or in other words row major indexing.
Because of that, it also makes sense to store the mat row-major, that is, laid out in memory as [1,2,3,4,5,6,7,8,9] so that mat[row][col] becomes mem[row*3 + col].
This is also how C does it. I think it was a mistake for OpenGL to break with that C convention given how close the relationship between C and OpenGL is. Mathematicians be damned.
>Some languages have fancy syntax that omits the need for braces. But for those that don't, this is clearly an array of 3 arrays of 3 ints.
>It's only intuitive to think of this as an array of rows given how it's laid out visually.
That's often the least important thing when designing a fast library. It's much better to optimize for common operations in your problem to be fast and thus choose the storage layout of a individual matrix instance accordingly (row-major, column-major, triangular, diagonal, block diagonal, band) so your data cache doesn't get tripped up regularily.
The actual API for creating new matrices can still be standardized and always the same. You can have the standardized API accept four column vectors for construction just as well. Just document it as such.
If you want, you can also represent the different kinds of vectors as different types, so you don't get identification problems that shouldn't exist. Also, lists of vectors are just that and there's no reason to represent both "list" and "vector" as "array". There are reasons to make the API for "vector" richer--especially since that's why the term "vector" was defined in the first place.
But for tiny 4x4 matrices, who cares about caching. But 4x4 are a minority use of matrices. Try 4096x4096.
As for OpenGL, they chose to keep the internal storage of matrices the same as it was on IRIS GL. That's a good decision for OpenGL and orthogonal to the question of how column vectors of that matrix can be accessed via an API and to the question of how row vectors of that matrix can be accessed via an API and to the question of how matrices can be created via an API.
From a mathematics standpoint, a matrix represents a linear map (a map is a function that takes a vector in space V to a vector in space W). It makes sense to be able to use matrix and linear map interchangeably as a user (i.e. use them just like you would call functions).
>It's only intuitive to think of this as an array of rows given how it's laid out visually.
That's often the least important thing when designing a fast library. It's much better to optimize for common operations in your problem to be fast and thus choose the storage layout of a individual matrix instance accordingly (row-major, column-major, triangular, diagonal, block diagonal, band) so your data cache doesn't get tripped up regularily.
The actual API for creating new matrices can still be standardized and always the same. You can have the standardized API accept four column vectors for construction just as well. Just document it as such.
If you want, you can also represent the different kinds of vectors as different types, so you don't get identification problems that shouldn't exist. Also, lists of vectors are just that and there's no reason to represent both "list" and "vector" as "array". There are reasons to make the API for "vector" richer--especially since that's why the term "vector" was defined in the first place.
But for tiny 4x4 matrices, who cares about caching. But 4x4 are a minority use of matrices. Try 4096x4096.
As for OpenGL, they chose to keep the internal storage of matrices the same as it was on IRIS GL. That's a good decision for OpenGL and orthogonal to the question of how column vectors of that matrix can be accessed via an API and to the question of how row vectors of that matrix can be accessed via an API and to the question of how matrices can be created via an API.
From a mathematics standpoint, a matrix represents a linear map (a map is a function that takes a vector in space V to a vector in space W). It makes sense to be able to use matrix and linear map interchangeably as a user (i.e. use them just like you would call functions).
> Try 4096x4096
Where should I, if we keep the context of OpenGL?
Apart from that, I consider this an implementation issue. My post was mostly about keeping a consistent "front end" semantics. However the compiler or library transforms the data does not pertain to it.
To be fair, I mentioned memory layout. If I had left that out, would we agree with eachother?
Where should I, if we keep the context of OpenGL?
Apart from that, I consider this an implementation issue. My post was mostly about keeping a consistent "front end" semantics. However the compiler or library transforms the data does not pertain to it.
To be fair, I mentioned memory layout. If I had left that out, would we agree with eachother?
(Column|Row) Vectors and matrices are abstractions, and having abstractions that are always the same API in all of programming (not just OpenGL) would be nice.
Of course making them again and again different and dumb is easier. So we got a lot of programs with weird arrays with no information on what matrix is what style. And thus readers get confused because half the information is missing in the source code! That makes things needlessly obtuse.
But let's say one would try to make the API for column vectors (resp row vectors resp matrices) always the same (OpenGL or not), then it's important to look at bigger matrices. Otherwise you wouldn't care about storage layout at all and that's gonna come back to bite you in the ass (several orders of magnitude slower than possible).
>But when you define a matrix in most programming languages you can only really write it like this:
Of course making them again and again different and dumb is easier. So we got a lot of programs with weird arrays with no information on what matrix is what style. And thus readers get confused because half the information is missing in the source code! That makes things needlessly obtuse.
But let's say one would try to make the API for column vectors (resp row vectors resp matrices) always the same (OpenGL or not), then it's important to look at bigger matrices. Otherwise you wouldn't care about storage layout at all and that's gonna come back to bite you in the ass (several orders of magnitude slower than possible).
>But when you define a matrix in most programming languages you can only really write it like this:
mat = [
[1,2,3],
[4,5,6],
[7,8,9]
]
I wouldn't repurpose arrays like that.Very true, and I think this is generally why C 'does it right'. But there are caveats.
Firstly, almost no-one actually creates their matrices as literals. Instead they are build dynamically an usually incrementally. Nevertheless, people still think back to these very basic examples when trying to understand things, so this probably should be standard.
Secondly, these decisions matter a lot for matrix vector multiplication. The 'C' way is very efficient for computing Ax which presumes vector x to be a Column vector. When you want to compute xA instead which presumes vector x to be a row vector, you want column-major ordering for efficiency. The same works for matrix-matrix multiplication AB where you want A to be row-major and B to be column-major.
To sum up. Because of our writing direction, and the natural way to write 2-d arrays, it makes sense to have row-major ordering. That immediately means it also makes sense to use column vectors, generally speaking. There are also exceptions if you happen to know if a matrix will often be at the 'left' or 'right' side of a multiplication.
Some more interesting notes that are less relevant:
In C the matrix you defined has an incompatible type with how most matrices are stored (most people in C store rows contiguously, then define a matrix as an array of pointers to rows).
Another common ways to use matrices is element-wise, where you refer to x[i][j]. Mathematicians expect i to be the row and j to be the column. Which matches C. Whilst Fortran does it the other way around.
Firstly, almost no-one actually creates their matrices as literals. Instead they are build dynamically an usually incrementally. Nevertheless, people still think back to these very basic examples when trying to understand things, so this probably should be standard.
Secondly, these decisions matter a lot for matrix vector multiplication. The 'C' way is very efficient for computing Ax which presumes vector x to be a Column vector. When you want to compute xA instead which presumes vector x to be a row vector, you want column-major ordering for efficiency. The same works for matrix-matrix multiplication AB where you want A to be row-major and B to be column-major.
To sum up. Because of our writing direction, and the natural way to write 2-d arrays, it makes sense to have row-major ordering. That immediately means it also makes sense to use column vectors, generally speaking. There are also exceptions if you happen to know if a matrix will often be at the 'left' or 'right' side of a multiplication.
Some more interesting notes that are less relevant:
In C the matrix you defined has an incompatible type with how most matrices are stored (most people in C store rows contiguously, then define a matrix as an array of pointers to rows).
int mat[3][3]; vs int *mat[3]
These are meaningfully different. In the first case getting element `mat[1][1]` is just 'base-pointer + 4'. Whereas in the second case, getting element `mat[1][1]` is 'read pointer at base-pointer + 1 and add 1 to the added pointer`.Another common ways to use matrices is element-wise, where you refer to x[i][j]. Mathematicians expect i to be the row and j to be the column. Which matches C. Whilst Fortran does it the other way around.
> The 'C' way is very efficient for computing Ax which presumes vector x to be a Column vector. When you want to compute xA instead which presumes vector x to be a row vector, you want column-major ordering for efficiency.
You have this the wrong way around. The C row-major ordering is better for A^Tx, the Fortran column-major format is better for Ax.
If you do Ax in row-major, you end up multiplying a row of A by x as your simd vector op. This then leaves you needing to horizontally reduce the result, which you avoid by using a column-major layout that instead accumulates multiple different results.
You have this the wrong way around. The C row-major ordering is better for A^Tx, the Fortran column-major format is better for Ax.
If you do Ax in row-major, you end up multiplying a row of A by x as your simd vector op. This then leaves you needing to horizontally reduce the result, which you avoid by using a column-major layout that instead accumulates multiple different results.
edit: My argument was based on the naive implementation of matrix-vector multiplication. But, as explained here [1]
there are more efficient algorithms for small-ish matrices that flip this around. For bigger matrices it remains better to do the efficient algorithm block-wise.
There is an argument that, for 'obviousness', the naive method matter more than the efficient one. But that gets tenuous. Especially because with some exceptions (finite fields, ugh) people are probably much better of using libraries for matrix operations, rather than writing their own.
[1] https://news.ycombinator.com/item?id=33359654
There is an argument that, for 'obviousness', the naive method matter more than the efficient one. But that gets tenuous. Especially because with some exceptions (finite fields, ugh) people are probably much better of using libraries for matrix operations, rather than writing their own.
[1] https://news.ycombinator.com/item?id=33359654
wait, so a row major stored matrix should be multiplied by extracting the columns first? Or what?
Is there some sort of animation that illustrates Ax and xA for both row major storage and column major storage?
Is there some sort of animation that illustrates Ax and xA for both row major storage and column major storage?
The schoolbook formulae for most linear algebra operations are written using scalar products, e.g. for Ax as scalar products of rows of A with x.
On a computer, such operations should almost never be implemented with scalar products, because the speed of a scalar product is limited by the latency of the FMA operation, not by its throughput. It is possible to accelerate the computation of a scalar product by ordering the elementary operations in a tree, but that complicates the program and it still does not reach the maximum speed of the hardware.
Except for vector-vector operations (i.e. for level 2 BLAS or greater levels), it is possible to change the order of the loops so that the scalar products are replaced by AXPY operations (whose result vectors should be kept in registers, to not be limited by the memory transfer throughput), which are not limited by the latency of the FMA, like the scalar products.
Except for vector-vector operations and matrix-vector operations (i.e. for level 3 BLAS or greater levels), it is possible to change the order of the loops so that the scalar products are replaced by rank-one matrix updates.
Most of the computation of a rank-one matrix update consists of the tensor product of 2 vectors.
Therefore, a matrix-matrix multiplication can be computed either by scalar products of the rows of the 1st matrix with the columns of the 2nd matrix, or by the tensor products of the columns of the 1st matrix with the rows of the 2nd matrix. The second method is much faster. The result of the tensor product must be kept in registers for the duration of the computation, so large matrices are partitioned in small blocks, i.e. submatrices, which are multiplied directly (an additional complication that increases the number of nested loops is that the small blocks that can be multiplied directly must be grouped in larger blocks that can be kept in the level-2 cache memory and reused for computations).
For matrix-matrix multiplications it is less important whether they are stored in column major order or in row major order, because when submatrices are brought into the L2 cache, the matrix elements will be gathered or scattered to be in the best order, before being loaded into registers repeatedly.
For matrix-vector multiplications, it is better to have the matrix in column major order, in order to compute the product by AXPY operations between columns and elements of the vector, and not by scalar products between rows and the vector (computing AXPY operations with shorter parts of a column at a time, so that the result vector can be kept in registers).
On a computer, such operations should almost never be implemented with scalar products, because the speed of a scalar product is limited by the latency of the FMA operation, not by its throughput. It is possible to accelerate the computation of a scalar product by ordering the elementary operations in a tree, but that complicates the program and it still does not reach the maximum speed of the hardware.
Except for vector-vector operations (i.e. for level 2 BLAS or greater levels), it is possible to change the order of the loops so that the scalar products are replaced by AXPY operations (whose result vectors should be kept in registers, to not be limited by the memory transfer throughput), which are not limited by the latency of the FMA, like the scalar products.
Except for vector-vector operations and matrix-vector operations (i.e. for level 3 BLAS or greater levels), it is possible to change the order of the loops so that the scalar products are replaced by rank-one matrix updates.
Most of the computation of a rank-one matrix update consists of the tensor product of 2 vectors.
Therefore, a matrix-matrix multiplication can be computed either by scalar products of the rows of the 1st matrix with the columns of the 2nd matrix, or by the tensor products of the columns of the 1st matrix with the rows of the 2nd matrix. The second method is much faster. The result of the tensor product must be kept in registers for the duration of the computation, so large matrices are partitioned in small blocks, i.e. submatrices, which are multiplied directly (an additional complication that increases the number of nested loops is that the small blocks that can be multiplied directly must be grouped in larger blocks that can be kept in the level-2 cache memory and reused for computations).
For matrix-matrix multiplications it is less important whether they are stored in column major order or in row major order, because when submatrices are brought into the L2 cache, the matrix elements will be gathered or scattered to be in the best order, before being loaded into registers repeatedly.
For matrix-vector multiplications, it is better to have the matrix in column major order, in order to compute the product by AXPY operations between columns and elements of the vector, and not by scalar products between rows and the vector (computing AXPY operations with shorter parts of a column at a time, so that the result vector can be kept in registers).
Thanks for this excellent explanation. I had questions, but they were answered by a closer reading of your explanation.
I recently did some work where, interestingly enough, this did not apply. I was doing linear algebra modulo 2 (i.e. binary entries, adding is xor, multiplication is and). In order to have any kind of memory (throughput) efficiency, you want to pack your matrix entries into bits of integers so a single 64 int can represent 64 entries.
This packing (when done row-wise or column-wise) gives a very quick primitive for computing inner products. <a, b> is just (in pseudo-code)
I wonder if a similar very fast scalar product primitive for floats (or ints?) could work. If not, I imagine that is due to FMA having some inherent efficiencies that cannot be gained by a very quick accumulation step. One of the beauties in the binary case is the popcnt instruction doing a 64-input accumulation very quickly.
I recently did some work where, interestingly enough, this did not apply. I was doing linear algebra modulo 2 (i.e. binary entries, adding is xor, multiplication is and). In order to have any kind of memory (throughput) efficiency, you want to pack your matrix entries into bits of integers so a single 64 int can represent 64 entries.
This packing (when done row-wise or column-wise) gives a very quick primitive for computing inner products. <a, b> is just (in pseudo-code)
sum(popcnt(a[i] & b[i]) for i in range len(a)) % 2
Hence in matrix multiplication over binary matrix you very much do want the l.h.s. to be row-major and the r.h.s. to be column major. And, in general, you very much want to read and write matrices in a contiguous manner. Because otherwise, depending on alignment, you are accessing 8, 32, or 64 bits of memory per single bit entry of matrix your are reading or writing.I wonder if a similar very fast scalar product primitive for floats (or ints?) could work. If not, I imagine that is due to FMA having some inherent efficiencies that cannot be gained by a very quick accumulation step. One of the beauties in the binary case is the popcnt instruction doing a 64-input accumulation very quickly.
> wait, so a row major stored matrix should be multiplied by extracting the columns first? Or what?
No, the other way around. In row-major storage, all the elements of the same row are contiguous, so it is faster to process the array row by row.
> Is there some sort of animation that illustrates Ax and xA for both row major storage and column major storage?
The explanations here are decent enough. All you need is the figure at the top, really: https://en.m.wikipedia.org/wiki/Row-_and_column-major_order .
No, the other way around. In row-major storage, all the elements of the same row are contiguous, so it is faster to process the array row by row.
> Is there some sort of animation that illustrates Ax and xA for both row major storage and column major storage?
The explanations here are decent enough. All you need is the figure at the top, really: https://en.m.wikipedia.org/wiki/Row-_and_column-major_order .
> Another common ways to use matrices is element-wise, where you refer to x[i][j]. Mathematicians expect i to be the row and j to be the column. Which matches C. Whilst Fortran does it the other way around.
That’s not quite it. The fact that the row index is first is purely a convention. And by convention, i in In X(i,j) in some Fortran code is still the row index.
The difference with C is the layout in memory, not the semantics: the fact that X(i+1,j) is contiguous to X(i,j) in Fortran and not in C, which means that for efficiency reasons inner loops in Fortran are over row indices rather than column indices as in C.
It would have been odd for Fortran, which has been since the beginning a language for numerical calculations, to “not match what mathematicians expect”.
That’s not quite it. The fact that the row index is first is purely a convention. And by convention, i in In X(i,j) in some Fortran code is still the row index.
The difference with C is the layout in memory, not the semantics: the fact that X(i+1,j) is contiguous to X(i,j) in Fortran and not in C, which means that for efficiency reasons inner loops in Fortran are over row indices rather than column indices as in C.
It would have been odd for Fortran, which has been since the beginning a language for numerical calculations, to “not match what mathematicians expect”.
vs
int (*mat)[3] = calloc(3, sizeof(int[3]));
which is a pointer to (or array of) arrays of 3 ints. Laid out flat, and you still get the `mat[i][j]` notation.That is indeed almost the same as
int mat[3][3]
The biggest reason people use the separate array of pointers is for very quick row swaps. Often comes up in matrix operations when trying to solve systems or put matrices in triangular forms.I think that whatever combination of tools/docs/compilers that makes row vs. column major notationally equivalent, amenable to introspection, and equally amenable to optimization will win big.
I’ve written godawful python extensions to feed TF or PyTorch mmapped regions to the pipeline in an efficient way. This is currently the way of things.
But half the trouble in model architecture engineering is keeping the ranks straight and having tools to debug when you messed up.
Einstein struggled with tensor math, read his letters. But he didn’t have modern computers or autodiff either. For that era, thinking in higher dimensions was a job for a genius.
I harbor a vague suspicion that some people get paid more because it’s needlessly hard.
I’ve written godawful python extensions to feed TF or PyTorch mmapped regions to the pipeline in an efficient way. This is currently the way of things.
But half the trouble in model architecture engineering is keeping the ranks straight and having tools to debug when you messed up.
Einstein struggled with tensor math, read his letters. But he didn’t have modern computers or autodiff either. For that era, thinking in higher dimensions was a job for a genius.
I harbor a vague suspicion that some people get paid more because it’s needlessly hard.
einops, einsum, pytrees, named tensors, and xla/jit/tensor compilation (jax/tvm/taco) will hopefully improve things
Hey! Don't just blame the mathematicians, if it were up to them you could notate column vectors as actual columns.
Mathematicians use row-major indexing just like any sensible person.
Mathematicians use row-major indexing just like any sensible person.
Mathematicians index the ith row and jth column of their n x m matrix as (i, j), which is nothing to do with row-major or column-major order. Those orderings are how a two-dimensional index is flattened to a one-dimensional index. Row-major makes it so that each row is stored one after the other: (i, j) goes to mi + j. Column-major uses i + nj instead.
That's as may be but M[i][j] and Mij still refer to the same element.
This doesn't tell us anything about memory layout though.
If a mathematician actually has to write out their matrices, they do it on a piece of paper. Because it is a 2D medium, the each element can have two "neighbors" on in the next row, and one in the next column.
If they do their math on a computer, they have the same problem as the rest of us -- but, they probably just use BLAS, which is typically column major. Although, I think some of the more modern BLAS spinoffs accept a layout option:
https://www.intel.com/content/www/us/en/develop/documentatio...
If a mathematician actually has to write out their matrices, they do it on a piece of paper. Because it is a 2D medium, the each element can have two "neighbors" on in the next row, and one in the next column.
If they do their math on a computer, they have the same problem as the rest of us -- but, they probably just use BLAS, which is typically column major. Although, I think some of the more modern BLAS spinoffs accept a layout option:
https://www.intel.com/content/www/us/en/develop/documentatio...
> Abstractly, a vector is neither a row nor a column. An n-dimensional vector is an ordered collection of n coordinates. An m by n matrix is an m-dimensional vector of n-dimensional vectors.
Oh we can get more abstract than that.
A vector is an element of a vector space, or, if you'd like, an F-module where F is a field. A matrix is a homomorphism between vector spaces.
Oh we can get more abstract than that.
A vector is an element of a vector space, or, if you'd like, an F-module where F is a field. A matrix is a homomorphism between vector spaces.
To give an explicit example that breaks the column-vector mold: the set of 2x2 matrices with real entries forms a vector space, where each vector is a 2x2 matrix.
Let's get even weirder and leave finite dimensions behind. The set of sequences of complex numbers z_n such that the infinite sum from i = 1 to infinity of |z_i|^2 converges.
Yeah having used Matlab extensively it's pretty clear that maths just got some things wrong. 1-based indexing is wrong. Matrix index ordering is wrong (try dealing with images in Matlab if you don't believe me).
Mathematicians generally don't care because equations are basically pseudocode and they can make up new syntax and hand-wave edge cases away whenever they want.
Mathematicians generally don't care because equations are basically pseudocode and they can make up new syntax and hand-wave edge cases away whenever they want.
It's common in mathematics to consider vectors as column vectors because the action of some linear operator represented by matrix A on a vector v is just the matrix product A*v.
Mathematical conventions are really dependent on each other. If you change one — a lot of others become awkward.
I wrote short notes about conventions in mathematics to understand them better: (pdf) https://github.com/dandanua/little-endian-vs-big-endian-in-q...
Mathematical conventions are really dependent on each other. If you change one — a lot of others become awkward.
I wrote short notes about conventions in mathematics to understand them better: (pdf) https://github.com/dandanua/little-endian-vs-big-endian-in-q...
The because seems wrong. I would argue vectors are columns because systems of linear equations are sentences which take up a row in a math paper.
the matrix is simply the layout and so R^n is column on the other side of the equals. It's arbitrary.
the matrix is simply the layout and so R^n is column on the other side of the equals. It's arbitrary.
Well, from the historical point of view I think you're right.
But in an abstract sense, matrix is a presentation of a linear operator in a specific basis. If we choose to write it as a transposed matrix, for example, then the matrix product rule would be different (it should correspond to the composition of linear operators).
But in an abstract sense, matrix is a presentation of a linear operator in a specific basis. If we choose to write it as a transposed matrix, for example, then the matrix product rule would be different (it should correspond to the composition of linear operators).
But I think this mixes issues.
Vectors are objects from some set (suppose R^n), linear operators some other set (suppose R^n -> R^n).
There exists an operation function application which takes an operator and a vector and produces a vector. There is an operation function composition which takes 2 operators to make a third.
There is a "functor" from vectors to coordinate vectors, linear operators to matrices, function application to multiplication, function composition to multiplication.
The rule for matrix multiplication is constrained by the functor. So you can't choose to write it incorrectly or its not a "homomorphism".
This is what I'm trying to express. I don't think you are "free to choose" in some sense. And I wouldn't say the orientation of the vector is caused by the multiplication rule over matrices. I think the "functor" is what is doing the choosing.
Vectors are objects from some set (suppose R^n), linear operators some other set (suppose R^n -> R^n).
There exists an operation function application which takes an operator and a vector and produces a vector. There is an operation function composition which takes 2 operators to make a third.
There is a "functor" from vectors to coordinate vectors, linear operators to matrices, function application to multiplication, function composition to multiplication.
The rule for matrix multiplication is constrained by the functor. So you can't choose to write it incorrectly or its not a "homomorphism".
This is what I'm trying to express. I don't think you are "free to choose" in some sense. And I wouldn't say the orientation of the vector is caused by the multiplication rule over matrices. I think the "functor" is what is doing the choosing.
Right, the matrix product rule comes from our choice of the representation of a linear operator (it's what you call a functor). But we are free to choose this functor. I didn't want to complicate things in a comment, the referenced notes have this moment.
Can anyone explain why in physics you use a 3x3 matrix for rotation but in computer graphics you have a 4x4 with “1” in the lower right corner of the matrix.
These are called homogeneous coordinates, and as expected since you have 4x4 matrices, vectors are also 4D, usually in the form (x,y,z,1). The idea is that you can multiply each coordinate by a non-zero scalar and the vector represents the same point, it has several interesting mathematical properties.
For computer graphics, the main advantage is that you can represent translation and projection in the matrix in addition to the other linear transformations like rotation and scaling.
For projection, the 4th vector coordinate, usually named "w", is called the homogeneous coordinate and it is typically not 1 after going through the projection matrix. As a final step, it divides all the other coordinates, so (x,y,z,w) becomes (x/w,y/w,z/w,1), x/w and y/w are the 2D screen coordinates, and z/w goes into the Z-buffer.
For computer graphics, the main advantage is that you can represent translation and projection in the matrix in addition to the other linear transformations like rotation and scaling.
For projection, the 4th vector coordinate, usually named "w", is called the homogeneous coordinate and it is typically not 1 after going through the projection matrix. As a final step, it divides all the other coordinates, so (x,y,z,w) becomes (x/w,y/w,z/w,1), x/w and y/w are the 2D screen coordinates, and z/w goes into the Z-buffer.
To add to the explanations already given here, the way I think of it as a non-maths person is:
A matrix represents a transform of the form:
Although this can represent scales and rotations, there's no way to represent a simple translation with this.
As proof, imagine you want to move everything 3 units along the x axis. You really just want to add 3 to x (x' = x + 3) but you can't: x' is always defined in terms of x, y and z (x' = Ax + By + Cz). There's no room for a constant.
To represent translations, then, what you really want is (x' = Ax + By + Cz + D), where D isn't multiplied by any component of the input vector, it's just D, your translation.
Well, it turns out you can do this by just adding an extra column to the matrix and using 1 for the fourth component of your vectors.
Now x' = Ax + By + Cz + Dw, where w=1, and D is your translation amount.
The full matrix then becomes
A matrix represents a transform of the form:
x' = Ax + By + Cz
y' = Dx + Ey + Fz
z' = Gx + Hy + Iz
...the A...I letters being the elements of the 3x3 matrix. (If you squint you can see the matrix above).Although this can represent scales and rotations, there's no way to represent a simple translation with this.
As proof, imagine you want to move everything 3 units along the x axis. You really just want to add 3 to x (x' = x + 3) but you can't: x' is always defined in terms of x, y and z (x' = Ax + By + Cz). There's no room for a constant.
To represent translations, then, what you really want is (x' = Ax + By + Cz + D), where D isn't multiplied by any component of the input vector, it's just D, your translation.
Well, it turns out you can do this by just adding an extra column to the matrix and using 1 for the fourth component of your vectors.
Now x' = Ax + By + Cz + Dw, where w=1, and D is your translation amount.
The full matrix then becomes
x' = Ax + By + Cz + Dw
y' = Ex + Fy + Gz + Hw
z' = Ix + Jy + Kz + Lw
w' = Mx + Ny + Oz + Pw
You can see how (D, H, L) now functions as a translation vector.Just to add to the existing responses, matrices represent linear transforms, which have many convenient properties. However, translation, that is moving things, is not a linear transform. (All linear transforms must map the zero vector to itself, but if you move an object at the origin somewhere else, it's no longer at the origin.) Translation is an affine transform.
Now here's the trick. You can embed your three-dimensional space into a four-dimensional space at a fixed fourth coordinate (but not zero). In this new space, translation in the original three-dimensional space is a linear transform.
Now here's the trick. You can embed your three-dimensional space into a four-dimensional space at a fixed fourth coordinate (but not zero). In this new space, translation in the original three-dimensional space is a linear transform.
The 4x4 represents both a 3x3 matrix A and a translation t, so that the “affine transformation” applied to a vector is v -> Av + t. One could pass around pairs (A, t) instead, and figure out how to compose those pairs, but it turns out there’s a nice way of embedding them in a 4x4 matrix such that composition is just 4x4 matrix multiplication.
You may opt into some other mathematical niceties too, which are handy for applying perspective transformations, but the main takeaway is that you get efficient compositions of affine transformations.
You may opt into some other mathematical niceties too, which are handy for applying perspective transformations, but the main takeaway is that you get efficient compositions of affine transformations.
Such 4x4 matrices allow for translations along with rotations (and skewing I think) into one linear transformation the 4th column can include the translation vector x,y,z and the last element is 1
[deleted]
> Recent mathematical treatments of linear algebra and related fields invariably treat vectors as columns (there are some technical reasons for this)
As someone who is studying the subject, I would be greatful if someone can please explain the technical reasons for this. I do not understand why I can't write them as row vectors and do the same Linear Algebra trickery (which we do with column vectors) another way around.
As someone who is studying the subject, I would be greatful if someone can please explain the technical reasons for this. I do not understand why I can't write them as row vectors and do the same Linear Algebra trickery (which we do with column vectors) another way around.
Matrix-vector multiplication works correctly when the vector is an Nx1 matrix, but if it's a 1xN matrix you need to take its transpose.
That's not just a mechanical computational consideration. One important interpretation of matrix-matrix multiplication is composition of linear transformations. When multiplying AB, you apply multiply matrix A by each column (vector) of B, and concatenate them to form the columns of the result AB. This corresponds to applying the transformation represented by A to a set of basis vectors. This correspondence is essential to the interpretation of matrix multiplication as a linear transformation, so it's important that it work out without additional manipulation.
Matrix-vector multiplication is more or less the heart of linear algebra, so in general it needs to be an elegant sensible operation.
That's not just a mechanical computational consideration. One important interpretation of matrix-matrix multiplication is composition of linear transformations. When multiplying AB, you apply multiply matrix A by each column (vector) of B, and concatenate them to form the columns of the result AB. This corresponds to applying the transformation represented by A to a set of basis vectors. This correspondence is essential to the interpretation of matrix multiplication as a linear transformation, so it's important that it work out without additional manipulation.
Matrix-vector multiplication is more or less the heart of linear algebra, so in general it needs to be an elegant sensible operation.
This is a little confused.
Matrix multiplication is defined as follows: if A is an m x n matrix and B is an n x p matrix, then the (i, j)th entry of C = AB is:
\sum_{k=1}^p A_{ik} B_{kj}
This definition works whether you think of a column vector as a vector and a row vector as a linear functional or vice versa.
How you interpret "what's happening" (whether you're expanding in a basis or taking dot products with basis vectors) depends on how you interpret the rows and columns of a matrix, regardless of the "orientation".
Matrix multiplication is defined as follows: if A is an m x n matrix and B is an n x p matrix, then the (i, j)th entry of C = AB is:
\sum_{k=1}^p A_{ik} B_{kj}
This definition works whether you think of a column vector as a vector and a row vector as a linear functional or vice versa.
How you interpret "what's happening" (whether you're expanding in a basis or taking dot products with basis vectors) depends on how you interpret the rows and columns of a matrix, regardless of the "orientation".
> depends on how you interpret the rows and columns of a matrix
I assume you know everything I am writing below, but maybe it's a clearer version of what I was trying to say in my other post.
Computationally, the definition of matrix multiplication AB = C is identical to multiplying A by each column of B, and concatenating the results as columns of C. The definition works both ways as you said, but this specific interpretation is one consequence of that definition.
Ax is a linear combination of the columns of A, with the values of x as coefficients. This is a direct translation of plugging the x values into a system of linear equations where the coefficients of each equation in the system form the rows of A:
You can go on like this, but the overall idea is that if you define "a matrix" such that each row contains the coefficients of one equation in a system of linear equations, then you are setting up row vectors to correspond to transformations (as in, y = a1x1 + a2x2), and then it becomes very ergonomic and natural to represent vectors as column vectors.
In short: if transformations go on the left and vectors go on the right, then transformations are rows and vectors are columns.
I assume you know everything I am writing below, but maybe it's a clearer version of what I was trying to say in my other post.
Computationally, the definition of matrix multiplication AB = C is identical to multiplying A by each column of B, and concatenating the results as columns of C. The definition works both ways as you said, but this specific interpretation is one consequence of that definition.
Ax is a linear combination of the columns of A, with the values of x as coefficients. This is a direct translation of plugging the x values into a system of linear equations where the coefficients of each equation in the system form the rows of A:
a11 * x1 + a12 * x2 = y1
a21 * x1 + a22 * x2 = y2
[ a11 a12 ; a21 a22 ] [ x1 ; x2 ] = [ y1 ; y2 ]
Moreover, in a matrix that corresponds to a system of non-redundant linear equations, the columns correspond to a basis for the image of the linear transformation represented by A. This follows immediately from the definition of a basis.You can go on like this, but the overall idea is that if you define "a matrix" such that each row contains the coefficients of one equation in a system of linear equations, then you are setting up row vectors to correspond to transformations (as in, y = a1x1 + a2x2), and then it becomes very ergonomic and natural to represent vectors as column vectors.
In short: if transformations go on the left and vectors go on the right, then transformations are rows and vectors are columns.
There’s no real reason to think of column vectors as points in space rather than the other way around when you’re doing linear algebra. It’s just a convention. If you wanted, you could think of a row vector as a point in space and transpose everything. In that case, your dual space of linear functionals would consist of column vectors rather than row vectors.
I second your request for a detailed explanation of the "technical reason" -- different people kept mentioning it, then saying that they would omit it. I want to know the technical / mathematical details for why, in isolation, row(1, 2, 3) is different from column(1, 2, 3). Like, all the way back down to the purest axioms or geometric intuitions.
I understand why the "rules" of matrix multiplication make the distinction -- but what is underneath those rules? What does it mean?
I understand why the "rules" of matrix multiplication make the distinction -- but what is underneath those rules? What does it mean?
Say you consider a column vector to be an actual physical vector in space. The geometric way to think of this entity is as a direction (n-1 degrees of freedom) and a magnitude (1 DOF).
Because of the definition of matrix multiplication, if we have two column vectors of length n, u and v, and u' denotes the transpose of u, then:
For this reason, it makes sense to think of column vectors as "vector vectors" and row vectors as linear transformations, which are not the same thing. Sometimes these linear transformations are called "linear functionals", but this point isn't particularly important. The important distinction is between "geometric objects in space" and "operators which encode actions on that space". You can think of these things as having different types.
To actually answer your question: in the above, there is no reason I couldn't have decided that row vectors are "vector vectors", resulting in column vectors being thought of as linear functionals. There is no technical reason motivating this distinction as far as I know. Column vectors are simply the convention in mathematics, while computer graphics did it the other way around for a while.
Because of the definition of matrix multiplication, if we have two column vectors of length n, u and v, and u' denotes the transpose of u, then:
u'v = |u| |v| cos \theta
where \theta is the angle between the directions determined by the vectors u and v. This means that the matrix u' (a row vector) can be thought of as a linear transformation mapping v to a measurement which encodes the angle between u and v (it's proportional to the scalar projection onto u).For this reason, it makes sense to think of column vectors as "vector vectors" and row vectors as linear transformations, which are not the same thing. Sometimes these linear transformations are called "linear functionals", but this point isn't particularly important. The important distinction is between "geometric objects in space" and "operators which encode actions on that space". You can think of these things as having different types.
To actually answer your question: in the above, there is no reason I couldn't have decided that row vectors are "vector vectors", resulting in column vectors being thought of as linear functionals. There is no technical reason motivating this distinction as far as I know. Column vectors are simply the convention in mathematics, while computer graphics did it the other way around for a while.
This is getting close to what I wanted -- thank you. But now you introduce what might be the real heart of it: the definition of matrix multiplication. I know the rules of this operation, and how matrix multiplication distinguishes rows and columns, and the geometric interpretation of a vector multiplied by a matrix. But where does that definition come from in the first place?
This might be too ill-posed to answer, it's been years since I was familiar enough with everything to know how to express exactly my confusion. But if you (or anyone) wants to take one more crack, I will appreciate it :)
This might be too ill-posed to answer, it's been years since I was familiar enough with everything to know how to express exactly my confusion. But if you (or anyone) wants to take one more crack, I will appreciate it :)
I don't think I'm able to give a real crack at the "where does it come from" question, but an example that might be intuitive -- if you had a system of equations, like high school algebra stuff:
a w + b v = y
c w + d v = z
If you know w and v, and want y and z, this is equivalent to a matrix multiplication: [a b] [v] = [y]
[c d] [w] [z]
And if you are solving that classic highschool problem "you know y and z, what are v and w," it is equivalent to the linear system solve operation, or multiplication by the inverse of the matrix. [v] = ([a b])^-1 [f]
[w] ([c d]) [g]
I apologize for the formatting. Imagine any time I've got parens or brackets on multiple lines, they are their big equivalent that envelopes both lines.Matrix multiplication corresponds to the composition of related linear operators. So, it goes down to the representation of a linear operator by a matrix. But we have freedom here. For example, we could choose to represent a linear operator by the transpose matrix (of the usual). But this requires changing the matrix product rule as well (to have the meaning of composition).
By convention, (roughly), row vectors are functions and column vectors are data.
[1, 2, 3][2
4
6]
Means 1*2 + 2*4 + 6*3Honestly?
I think it's because column vectors are less of a pain to write on the blackboard.
I'm being serious. There is no reason to say that columns are vectors and rows are linear functionals. The other way around is a perfectly valid convention. But I did find that when I taught linear algebra, column vectors were just easier to write on the board.
I think it's because column vectors are less of a pain to write on the blackboard.
I'm being serious. There is no reason to say that columns are vectors and rows are linear functionals. The other way around is a perfectly valid convention. But I did find that when I taught linear algebra, column vectors were just easier to write on the board.
You can, but you'll be doing linear algebra in the dual space, with covectors, instead of in the regular space with vectors.
Yes and yes, to the Homogeneous coordinate being the first column or row.
That way, you can have sparse vectors, or sparse matrices. It doesn't matter how many dimensions are specified for your object, stored in memory. It gets transformed the same.
Have a point at the origin? [1] Done. It doesn't matter if it's a 2D, 3D, 7D.
Have a 0-vector? [0] Done.
Have a 2D point? [1, x, y]
Have a 3D point? [1, x, y, z]
Have a 1D vector? [0, x]
I honestly don't know why we don't all do it this way.
That way, you can have sparse vectors, or sparse matrices. It doesn't matter how many dimensions are specified for your object, stored in memory. It gets transformed the same.
Have a point at the origin? [1] Done. It doesn't matter if it's a 2D, 3D, 7D.
Have a 0-vector? [0] Done.
Have a 2D point? [1, x, y]
Have a 3D point? [1, x, y, z]
Have a 1D vector? [0, x]
I honestly don't know why we don't all do it this way.
If f: X -> Y and g: Y -> Z are functions, their composite is the function from X to Z which most mathematicians would write as "g o f" (typographically that should be a small centered circle, not the letter "o" as here), defined by
Anyway, consider how this works if f and g are linear transformations and X, Y, and Z are finite-dimensional vector spaces, as in linear algebra. So we can take X = F^p, Y = F^n, and Z = F^m, where F is the field of scalars, and vectors are just "arrays of numbers".
Or you can bite the bullet, write arguments to functions on the left, and then defaulting to row vectors would be okay. It's just notation.
(g o f)(x) = g(f(x)).
Notice that the function f, which is applied first, appears on the right and g, which is applied second, appears on the left. But if we draw a diagram of the composite, it looks like this: f g
X ----> Y ----> Z
In the diagram, the first applied function is on the left and the second applied function is on the right. Some people have felt that, for this reason, it would be better to write composites with the arguments on the left, like this: (x)f and [x(f)]g . However, most mathematicians follow the right-to-left convention (that's the way it seems to be in most programming languages, right?).Anyway, consider how this works if f and g are linear transformations and X, Y, and Z are finite-dimensional vector spaces, as in linear algebra. So we can take X = F^p, Y = F^n, and Z = F^m, where F is the field of scalars, and vectors are just "arrays of numbers".
f: F^p -> F^n, g: F^n -> F^m
You can represent f and g using matrices: f is represented by an n x p matrix B and g is represented by an m x n matrix A. In that case, function application is done by matrix multiplication. The composite g o f is represented by the product A B, which is an m x p matrix. Notice that the matrices are in the same order in "A B" as their functions are in "g o f". Inputs to g o f are p-dimensional vectors; where should a p-dimensional vector go in order that multiplication by A B should represent application of the function g o f? A B [p x 1 column vector] ---> okay
[1 x p row vector] A B ---> wrong: you can't multiply (1 x p) * (m x p)
In other words, once you commit to the standard convention for writing composites, you're committed to using column vectors as inputs to linear transformations written as matrices, or somewhere you'll probably have an inconsistency in the order of application.Or you can bite the bullet, write arguments to functions on the left, and then defaulting to row vectors would be okay. It's just notation.
As someone who's in love with OpenGL, man what a great piece and interesting piece of history.
Not in love with GL as someone that's written GL drivers and re-implemented GL on top of modern APIs, but....
for me, the GL style (whatever you decide to call it) seems to make sense from a performance POV. I often want one axis or the translation. It's a single SIMD load since vs if they were stored the other way I'd have to load every 4th float.
Also, while you do type them out in code different from math
for me, the GL style (whatever you decide to call it) seems to make sense from a performance POV. I often want one axis or the translation. It's a single SIMD load since vs if they were stored the other way I'd have to load every 4th float.
Also, while you do type them out in code different from math
x x x x
y y y y
z z z z
t t t t
If you squint you can view them as columns xaxis = [
.
.
.
.
]
yaxis = [
.
.
.
.
]
zaxis = [
.
.
.
.
]
translation = [
.
.
.
.
]
matrix = [xaxis, yaxis, zaxis, translation] // 4 columns side by sideAgree, this email format is also so easy to parse. It's like I'm dropped right into the action.
As someone else said, when using einops, einsum, named tensors (named dimensions), it doesn't really matter.
For most software frameworks, there is anyway not a diff between column vector vs row vector, it's just shape [dim] in any case, but it becomes relevant how you do matmul etc. with some other matrix [dim1, dim2] or some other tensor [dim1, dim2, dim3].
For most software frameworks, there is anyway not a diff between column vector vs row vector, it's just shape [dim] in any case, but it becomes relevant how you do matmul etc. with some other matrix [dim1, dim2] or some other tensor [dim1, dim2, dim3].
Fix the bug:
Is the matrix row major or column major?
https://bugfix-66.com/e74dcfb6a10cb71404f495922533d947e39fd7...
Is the matrix row major or column major?
https://bugfix-66.com/e74dcfb6a10cb71404f495922533d947e39fd7...
Interesting that this mentions people using row vectors in physics. In my physics classes (and math classes) vectors were always represented as columns.
For me a vector is just abstractly three or four numbers together. (Maybe with some defined transformation properties if I'm doing physics.) Whether a matrix is stored row-major or column-major is an implementation detail I try to hide.
In fact I used to make indexing errors all the time. But now when I need to deal with vectors and matrices, I look up or decide on the "correct" way, hide everything behind Vector and Matrix classes, implement dot and cross product once, and never fall back to the raw vector elements if possible.