Sometimes programing can be fun. Sometimes it isn’t. I particularly hate debugging.
Today I was programming a function that would compute the commutator of two square matrices in c++, which do not have a definite size in the program. I ended up implementing it in the following way, by using the pointers to the beginning of the matrices and the size of the matrix n.
void commutator(int n , double *a1, double *a2, double *b1)
{
for(int i=0; i<n; i++){
for(int j=0; j<n; j++){
*(b1+i+n*j)=0;
for(int k=0; k<n;k++){
*(b1+i+n*j)+= *(a2+i+n*k)*(*(a1+k+j*n))-*(a1+i+n*k)*(*(a2+k+j*n)) ;
}
}
}
}
The two matrices are a1, a2, and the result is placed on b1.
The placing of parenthesis is very crucial. This was the best solution I found after trying all kinds of other ways that were giving me trouble with the compiler. This is probably good for numerics, but lousy for preventing errors in the handling if one does not pay attention. What is fascinating to me is that the symbol * does two things: both multiplication and looking up the address of a pointer. And these can appear consecutive to one another.
Oh well. Just thought I’d share my messy code.
Why don’t you use a well-developed library to handle your matrix operations, e.g. http://exmat.sourceforge.net/ ? There is no point in reinventing the wheel.
Hi Matt:
Oh that’s easy to explain: I’m relearning how to program properly and I’m doing this as an exercise. Part of the reason I want to do some of these things myself has to do with optimization on multiprocessors and other tricks to get the computer to work more efficiently.
Dear David, what I like best about this posting is the CSS design of the embedded code – a part of WordPress, right?
If you thought that the readers would be able to learn what the star can actually mean in C plus plus, well, you have failed at least in my case. I only understand * as multiplication.
There exists a two-step way to define the commutator:
http://www.wolfram.com/products/mathematicahomeedition/?greetingsfrom=lubos-motl-reference-frame
commutator[a_, b_] := a.b – b.a;
The first step listed above is to buy Mathematica 7.0.1.
Oh, I see. Now I understand the pointer. The matrix is simply squeezed somewhere in the memory as a one-dimensional array, a sequence of mud.
So if you want to move by j rows and i columns, it’s like moving by i+n*j places. Also, b1 stands for the physical address of some memory point where the b1 array-matrix-mud-junk begins
And *(something) is the memory cell corresponding to the address calculated as “something”. So
*(b1+i+n*j)=0;
is the command that annihilates the content of a number at i-th column and j-th row of the array-matrix-mud-junk called b1 because the content of the parenthesis is the right address.
The size of the cell is “double” i.e. 32 bits or, more likely, 64 bits?
Funny. Thanks for your C++ lecture. It’s clearly somewhat closer to the machine code than the higher languages. I would love assemblers 25 years ago but it seems nicer to leave this machine-code-level work to others.
The code looks almost like C rather than C++, except for where the variables are declared. C++ is a significantly higher-level language than C, but it keeps the C parts around anyway. Most C++ code usually avoids pointers as much as possible, or at least wraps them in safer data structures. Apparently David likes to muck around with them. Usually the efficiency gain isn’t worth the effort, but I’ve never had to optimize parallel code, so probably he has a good reason.
Mathematica is painfully slow.
Hi Lubos:
If Mathematica weren’t so slow, it could be useful for moderate numerical simulations. As it is, solving differential equations with matrices is too slow to be useful and one needs to resort to something that more directly addresses the machine.
Oh by the way, C and C++ is notoriously famous for being able to produce code that is almost incomprehensible to people. Unlike postcript which is essentially always incomprehensible to people
Dear David,
you could often be surprised that Mathematica can solve many differential equations either analytically, or self-optimize the numerical methods to solve them and outclever your C++ code. I am not saying that this will always occur but sometimes it does.
Decades ago, I would also think that addressing the transistors of a machine directly would improve my efficiency in dealing with computers – essentially because anything else that is added upon the machine code is garbage created by inferior people.
I no longer think so. There’s a lot of useful structure in the higher languages and you are often rewarded by great gifts for this apparent slowness. And by the way, the machine code was invented by mortal humans, too.
Cheers
LM
Premature optimization is the root of all evil!
You’d probably find your own code easier to read and understand if you didn’t use both meanings for the *
There’s really not a good reason to use:
When you can write:
This de-messifies it a bit, especially when you have multiple array references.
Note that any modern optimizing compiler is going to generate the same code for pair of accesses (if you don’t believe it, generate assembly for both sets of code and compare).
Yuk. C++ I know not why people cause themself such pain. I’ll file it under ‘programming fetishes’
Haelfix…. “Premature optimization is the root of all evil!” Truer words were never spoken
Not sure if the optimizer does this already, but rather than multiplying and adding, you might be better off creating dummy pointers that you increment by n and/or 1.
rather than:
*(a2+i+n*k)
//prep it
ptr2 = a2+i;
//do bare minimum in loop
for loop {
ptr1 += 1;
ptr2 +=n;
*res += (*ptr1)*(*ptr2);
}
So there is no chance any superfluous multiplication is done.
You should also const the input matrix pointers in the prototype. This tells the reader that they are input.
Hi everyone:
Thanks for the pointers
. I’ll implement them in the code better.