Hello! On Tue, Feb 20, 2007 at 04:16:08PM -0300, Charlls Quarra wrote:
I know i'm new to the list and probably you feel you have to be so overly cynical and aggresive to make some sort of social standpoint.
I'm sorry, no offence was intended... My point was that GiNaCs design makes some operations inefficient (and/or ugly), so GiNaC is not proper tool to perform them. But apparently such operations are not necessary for your task. So I wonder what is exactly aggressive and/or cynical in my posts...
Oh well, first you said that it wasn't possible to do it in ginac and i should look somewhere else, now you speak about it being inefficient.
Sure it is possible _in principle_, as in "the determined Real Programmer can write FORTRAN programs in any language". As for me "being unable to do something in a natural and efficient way" and "impossible" is essentially the same thing.
I happen to have read that chapter that you RTFM me before asking about write access, and precisely on that page it says that copy is done lazily, so in fact the copy operation:
matrix sub = ex_to<matrix>( mat(1,2) );
should be "light-weight", at least if the matrix is really just a handle like the ex.
Matrices are unshareable by default (because they provide write access to data): /** Very common ctor. Initializes to r x c-dimensional zero-matrix. * * @param r number of rows * @param c number of cols */ matrix::matrix(unsigned r, unsigned c) : inherited(TINFO_matrix), row(r), col(c), m(r*c, _ex0) { setflag(status_flags::not_shareable); } so the manual is somewhat inaccurate...
now, when i do sub(3,4) = 15;
this is arguably doing a copy (again, according to the RTFMed chapter) because both sub and mat(1,2) point to the same object,
They don't. For illustration consider the following code: // matmat.cpp #include <iostream> #include <string> #include <ginac/ginac.h> #include <stdexcept> using namespace std; using namespace GiNaC; void test_ex() { symbol x("x"), y("y"); ex foo = pow(x,2)*y + pow(y,2)*x + pow(y,3) + pow(x,3); cout << foo << endl; ex bar = foo; bar += pow(x,2)*pow(y,2); cout << bar << endl; } void test_matr() { matrix A = ex_to<matrix>(symbolic_matrix(2, 2, string("a"))); cout << A << endl; matrix B = A; B(0, 0) = ex(1); cout << B << endl; } int main(int argc, char** argv) { test_ex(); test_matr(); return 0; } Let's run it in the debugger: $ cat matmat.gdb break matmat.cpp:12 break matmat.cpp:14 break matmat.cpp:15 break matmat.cpp:21 break matmat.cpp:23 break matmat.cpp:24 commands 1 print foo.bp.p.refcount end commands 2 print foo.bp.p.refcount print bar.bp.p.refcount end commands 3 print foo.bp.p.refcount print bar.bp.p.refcount end commands 4 print A.refcount end commands 5 print A.refcount print B.refcount print A.m print B.m end commands 6 print A.refcount print B.refcount print A.m print B.m end $ echo $CXXFLAGS -Wall -O2 -g -pipe -funroll-loops -ffast-math -fomit-frame-pointer -finline-limit=1200 -m32 -march=pentium4 -mmmx -msse2 -mfpmath=sse $ g++-4.1 $CXXFLAGS -o matmat matmat.cpp -lginac $ gdb --silent ./matmat Using host libthread_db library "/lib/tls/i686/cmov/libthread_db.so.1". (gdb) source matmat.gdb Breakpoint 1 at 0x80491fb: file matmat.cpp, line 12. Breakpoint 2 at 0x8049225: file matmat.cpp, line 14. Breakpoint 3 at 0x804932a: file matmat.cpp, line 15. Breakpoint 4 at 0x8049a25: file matmat.cpp, line 21. Breakpoint 5 at 0x8049c47: file matmat.cpp, line 23. Breakpoint 6 at 0x8049ca7: file matmat.cpp, line 24. (gdb) run Starting program: /afs/diastp.jinr.ru/user/varg/tmp/matmat Breakpoint 1, test_ex () at matmat.cpp:12 12 cout << foo << endl; $1 = 1 (gdb) cont Continuing. x^3+y*x^2+y^3+y^2*x Breakpoint 2, test_ex () at matmat.cpp:14 14 bar += pow(x,2)*pow(y,2); $2 = 2 $3 = 2 // there is only one copy with reference count 2 (gdb) Continuing. Breakpoint 3, test_ex () at matmat.cpp:15 15 cout << bar << endl; $4 = 1 $5 = 1 // now there are two different expressions (gdb) Continuing. x^3+y*x^2+y^3+y^2*x^2+y^2*x Breakpoint 4, test_matr () at matmat.cpp:21 21 cout << A << endl; $6 = 0 // matrices are not sharable, so reference counting is pointless (gdb) Continuing. [[a00,a01],[a10,a11]] Breakpoint 5, test_matr () at matmat.cpp:23 23 B(0, 0) = ex(1); $7 = 0 $8 = 0 $9 = {<std::_Vector_base<GiNaC::ex,std::allocator<GiNaC::ex> >> = { _M_impl = {<std::allocator<GiNaC::ex>> = {<__gnu_cxx::new_allocator<GiNaC::ex>> = {<No data fields>}, <No data fields>}, _M_start = 0x805d710, _M_finish = 0x805d720, _M_end_of_storage = 0x805d720}}, <No data fields>} $10 = {<std::_Vector_base<GiNaC::ex,std::allocator<GiNaC::ex> >> = { _M_impl = {<std::allocator<GiNaC::ex>> = {<__gnu_cxx::new_allocator<GiNaC::ex>> = {<No data fields>}, <No data fields>}, _M_start = 0x805a9f0, _M_finish = 0x805aa00, _M_end_of_storage = 0x805aa00}}, <No data fields>} // Adresses of actual storage (_M_start, _M_finish) are different for // A and B/ Both A and B are _different_ objects (gdb) Continuing. Breakpoint 6, test_matr () at matmat.cpp:24 24 cout << B << endl; $11 = 0 $12 = 0 $13 = {<std::_Vector_base<GiNaC::ex,std::allocator<GiNaC::ex> >> = { _M_impl = {<std::allocator<GiNaC::ex>> = {<__gnu_cxx::new_allocator<GiNaC::ex>> = {<No data fields>}, <No data fields>}, _M_start = 0x805d710, _M_finish = 0x805d720, _M_end_of_storage = 0x805d720}}, <No data fields>} $14 = {<std::_Vector_base<GiNaC::ex,std::allocator<GiNaC::ex> >> = { _M_impl = {<std::allocator<GiNaC::ex>> = {<__gnu_cxx::new_allocator<GiNaC::ex>> = {<No data fields>}, <No data fields>}, _M_start = 0x805a9f0, _M_finish = 0x805aa00, _M_end_of_storage = 0x805aa00}}, <No data fields>} // The address of storage did not change when we assigned something to B // (_M_start and _M_finish are the same is in the previous breakpoint) (gdb) Continuing. [[1,a01],[a10,a11]] Program exited normally. (gdb)
but the write operation on a handle triggers the copy. However, now that we are talking about _efficiency_, this is not unavoidable, at least from what refcounting refers. If one would do:
matrix sub = ex_to<matrix>( mat( 1 , 2 ) ); mat( 1 , 2 ) = 0; // <-- this kills the first reference sub(3, 4) = 15; // no copy should be involved anymore, there is a single reference to the underlying matrix mat( 1 , 2 ) = sub;
so in this case, Apparently no deep copy should be involved.
Unfortunately, you are wrong: for matrices there is a copy operation in every assignment. Usually this is not a big deal: the actual elements of matrices are ex (which are usually shareable), so copying is not that expansive (but still there are nxm handles to copy). Storing non-shareable objects (like matrices) in matrix is _very_ inefficient.
I actually wasnt trying to be fixed on doing this with matrix or ginac objects only other than esthetical/completeness reasons: I even could come around to consider doing it with std::vector after all (or even better, blitz), but that was off the point by several miles. I was really just trying to learn ginac by doing something that i would need commonly.
And I was trying to explain that one don't have to treat every object as an expression (like in CASes). If it is natural to store (and manipulate) your data as an array of expressions (or whatever else) -- just use that. There is no need to pack <your data type> it into yet another expression.
I was also intrigued if others don't seem to require multi indexing (i use it all the time on maple using lists of lists, but it doesnt seem to be efficient with deep stacking).
I think because STL containers provide this functionality. Best regards, Alexei -- All science is either physics or stamp collecting.