I agree with your suggestion of adding a FAQ about this issue, I was bitten by it the first time I used blitz, and quite a few people have asked about it in the list over time. Having said that, there is a very good explanation of how to implement matrix products, and more generally tensor operations using blitz in the manual, see section 3.11 Tensor notation. I think that approach provides a much more intuitive mechanism to implement the matrix product than the overloaded functions, and you only have to write it once, as it is already type independent. In summary, blitz provides a mechanism to write the matrix product as you do in the math definition, if A and B are matrices and you put the product in C the mathematical definition of the product is:
Cij = sum(Aik*Bkj)
where the sum is over k, the example code in blitz is (from the manual):
...
Array<float,2> A,B,C;
C = sum(A(tensor::i,tensor::k)*B(tensor::k,tensor::j),tensor::k);
...
The tensor object is defined inside the blitz namespace, and the above works for any type, not just float, so long as the product of two elements of that type is defined. Note that you can use the tensor objects as you would indices to write most tensor expressions just as you would on paper! Isn't blitz great?
Hope this helps.
Salut
Carlos A. Rega, PhD
Development Scientist
Malvern Instruments Ltd
Grovewood Rd, Malvern
WR14 1XZ
Tel: +44 (0)1684 892 456
Fax: +44 (0)1684 892 789
e-mail: carlos.rega@malvern.co.uk
"Marc Vinyes" <mvinyes@iua.upf.es> Sent by: blitz-support-bounces@oonumerics.org
07/05/2005 16:53
Please respond to mvinyes; Please respond to Support list for Blitz++
To: blitz-support@oonumerics.org
cc:
Subject: RE: [Blitz-support] matrix product (matrix dimensions posted againwithin the list)
Thank you very much Derrick Bass and specially Julian C. for you clear and
detailed explanation.
Now I see that I totally misunderstood the meaning of the * operator and I
thought it was the matrix "standard multiplication" (from a mathematical
point of view).
I suggest adding an entry to the FAQ because I think that this issue might
confuse other people in the future...
> Blitz Arrays are intended to perform elementwise math operations and that
> is
> what the operator*() does. You get a 2x1 result Array because that is
> what
> you allocated.
Finally, as my algorithms use the matrix standard multiplication quite
often and I appreciate them to be overloaded (to code fast), I have
decided to use octave's classes or newmat's library.
However I will certainly store your function for furher uses of blitz in
other problems where the standard-*-matrix overloading is less useful or
in a final efficient implementation of the algorithm.
Again, thank you.
Sincerely,
MarC
>The result Array A is filled by marching through Arrays B
> and C in elementwise order and performing multiplication. These are
> row-major Arrays, so we get 1x1=1 and 2x2=4 for the results that are
> stored
> into A. Of course, this is not the semantics you are looking for. What
> you
> want is to loop over the rows of B and the columns of C and perform inner
> products to get the resulting element values of A. You can do this with
> the
> proper use of Array slicing and summation.
>
> void matmatProduct(Array<int,2> A, Array<int,2> B, Array<int,2> C) {
> int rows = B.rows();
> int columns = C.columns();
> assert(B.columns() == C.rows()); // sanity check
> assert(A.rows()==rows && A.columns()==columns); // could reshape A here
> instead
> for (int i=0; i<rows; ++i)
> for (int j=0; j<columns; ++j)
> A(i,j) = sum(B(i,Range::all())*C(Range::all(),j));
> }
>
> Sorry there is no simple way to do this by overloading the standard *
> notation, but this function will give you the result you expected.
>
> Regards, Julian C.
>
> Dr. Julian C. Cummings Office: PB-111
> Caltech/CACR, MC 158-79 Phone: 626-395-2543
> 1200 E. California Blvd. Fax: 626-584-5917
> Pasadena, CA 91125
_______________________________________________
Blitz-support mailing list
Blitz-support@oonumerics.org
http://www.oonumerics.org/mailman/listinfo.cgi/blitz-support
______________________________________________________________________
This email and any files transmitted with it are confidential and
intended solely for the use of the individual or entity to whom
they are addressed.
If you have received this email in error please notify the
originator of the message.
Any views expressed in this message are those of the individual
sender, except where the sender specifies and with authority,
states them to be the views of Malvern Instruments.
______________________________________________________________________