|
|
Blitz Bugs : |
From: Julian C. Cummings (cummings_at_[hidden])
Date: 2004-03-15 21:40:26
Hi Sandy,
Blitz Arrays support elementwise math operations. Thus, C=A*B will loop
over the elements of Arrays A and B, multiplying the element values from A
and B and assigning the result to the corresponding element in Array C. It
is assumed that all three Arrays have the same shape (same extent in each
dimension). In your example, Array B is not the same shape as Array A.
Blitz acted as if it were and padded the missing elements with zeroes.
What you want is a true Matrix class that does standard matrix
multiplication. There is a Matrix class in blitz, but the binary operator*
defined for multiplying two Matrix objects does an elementwise
multiplication similar to the Array class. Perhaps this should be changed
or perhaps a special product() function should be written to do the standard
matrix-matrix multiplication. Suggestions?
There is also a TinyMatrix class in which the number of rows and number of
columns are template parameters and hence compile-time constants. There is
already a product() function defined for doing matrix-matrix multiplication
with this object. So you can do
TinyMatrix<double,4,4> A;
TinyMatrix<double,4,2> B,C;
C = product(A,B);
If you are doing a lot of matrix operations with small matrices, you might
want to look at the tvmet package (tvmet.sourceforge.net), which has a more
complete implementation of these tiny Vector and Matrix classes. It is in
the works to replace the blitz TinyVector and TinyMatrix classes with the
classes from tvmet in the future. The blitz Vector and Matrix classes
(non-Tiny) have sizes that are set at run time and are meant for larger
containers. These classes are somewhat incomplete, as they were added later
on in the blitz development and never really got mature. However, there's
been a lot of interest in having real Vector and Matrix classes available,
along with some standard linear algebra routines such as LU decomposition.
I hope to see some work done on this also. There are other C++ packages out
there focusing on this area, such as uBlas, MTL and TNT.
Regards, Julian C.
Dr. Julian C. Cummings
Staff Scientist, CACR/Caltech
(626) 395-2543
cummings_at_[hidden]
> -----Original Message-----
> From: blitz-bugs-bounces_at_[hidden]
> [mailto:blitz-bugs-bounces_at_[hidden]] On Behalf Of Sandy Jackson
> Sent: Sunday, March 14, 2004 9:25 PM
> To: blitz-bugs_at_[hidden]
> Subject: [Blitz-bugs] Matrix multiplication
>
>
> I am not sure what the what the folowing expression is
> suppose to return for type Array<double, 2>:
>
> C = A*B;
>
> If it is standard matrix multiplication, then, I believe that
> I have found a bug in the implementation. If it is not
> intednded to be standard matrix multiplication, what does it return ?
>
> The following short example illustrates the problem:
>
> // Inclusion of C++ standard library header files
>
> #include <iomanip>
> #include <iostream>
> #include <istream>
> #include <ostream>
> #include <fstream>
> using namespace std;
>
> // Inclusion of C++ non-standard library header files
> #include <blitz/array.h>
> using namespace blitz;
> typedef Array<double, 2> Matrix;
>
> int main()
> {
>
>
> Matrix A(4,4);
> A = 1,2,3,4,2,2,2,2,3,3,3,3,4,4,4,4;
>
> Matrix B(4,2);
> B = 1,2,1,2,1,2,1,2;
> Matrix C(4,4), Bo(4,4);
>
> C = 0;
> Bo = 0;
>
> Bo = A*B;
>
> for (int i = 0; i < A.rows(); i++ )
> for ( int j = 0; j < B.cols(); j++ )
> for( int k = 0; k < A.cols(); k++ )
> C(i,j) += A(i,k) * B(k,j);
>
> cout << "A = " << A << endl << " B = " << B << endl;
> cout << " My A*B = " << C << endl;
> cout << " Blitz A*B = " << Bo << endl;
>
> return 0;
> }
>
> I compiled the above code using gcc-3.3.1-3 for i-686-pc-cygwin.
>
> I get the following results.
>
> A = 4 x 4
> [ 1 2 3 4
> 2 2 2 2
> 3 3 3 3
> 4 4 4 4 ]
>
> B = 4 x 2
> [ 1 2
> 1 2
> 1 2
> 1 2 ]
>
> My A*B = 4 x 4
> [ 10 20 0 0
> 8 16 0 0
> 12 24 0 0
> 16 32 0 0 ]
>
> Blitz A*B = 4 x 4
> [ 1 4 3 8
> 2 4 2 4
> 1.5915e-312 6.10782e-260 6.61505e-314 8.7214e-312
> 0 0 0 0 ]
>
>
>
>
>
> =====
> May God Bless You
>
> Sandy R. Jackson
>
> "Because He lives, I can face tomorrow; Because he lives, all
> fear is gone. Because, I know He holds the future. And life
> is worth the living - just because He lives."
>
> William Gather, 1971 _______________________________________________
> Blitz-bugs mailing list
> Blitz-bugs_at_[hidden]
> http://www.oonumerics.org/mailman/listinfo.cgi> /blitz-bugs
>