|
|
Ublas : |
Subject: Re: [ublas] Does uBLAS work with dd_real and qd_real ?
From: Nasos Iliopoulos (nasos_i_at_[hidden])
Date: 2010-03-25 17:23:08
Paul,
you have some very interesting points here. Let me express my view also.
I don't think that all types that can be possibly represented in ublas include zero as an element Zero is more or less an arbitrary choice (which is fine if it is in a standard), but I don't think it should be adopted in a generic library. The underlying types' default constructor should be responsible for the default initialization. uBlas can only define a standard behaviour, and the less number of assumptions it makes the better.
Furthermore, since the standard requires POD number types (double, float ints etc.) to be zero-initialized when default constructed, why shouldn't this be adopted by other number types that include zero in their set definition?
My only scepticism with intialization in default construction relates to performance issues when creating large arrays or many small ones. Of course there is a simple work-around: create a custom-allocator that does not execute default initialization for POD types when you want construction performance.
It would be interesting to hear what other people think.
Best
Nasos
> From: paul.leopardi_at_[hidden]
> To: ublas_at_[hidden]
> Date: Fri, 26 Mar 2010 01:33:16 +1100
> Subject: Re: [ublas] Does uBLAS work with dd_real and qd_real ?
>
> On Friday 26 March 2010 00:27:09 Nasos Iliopoulos wrote:
> > Paul,
> > the problem seems to be in the implementation of the dd_real in the qd
> > library. I noticed that the default constructor does not initialize the
> > number class, hence constructing it in an ill-formed state. If you
> > recompile the qd library changing dd_real() { }
> > to:
> > dd_real() { x[0]=0.0; x[1]=0.0;}
>
> Hi Nasos,
> Yes, I noticed this too. To my mind this is strange behaviour on the part of
> uBLAS. Specifically, if the type of m is matrix<scalar_t> then m.clear() does
> not set the values of m to scalar_t(0) but rather to scalar_t(), i.e. it calls
> the default constructor rather than the constructor from 0 (if any).
>
> What does the C++ standard say about the result of calling the default
> constructor for float or double [ i.e. float(), double() ]? It says that since
> float and double are POD (Plain Old Data) classes, the result is zero-
> initialization. But for non-POD classes T, calling T() just results in default
> initialization.
>
> http://www.fnal.gov/docs/working-groups/fpcltf/Pkg/ISOcxx/doc/POD.html
> http://www.boost.org/doc/libs/1_33_1/libs/utility/value_init.htm
>
> Should dd_real be changed so that dd_real() yields 0, thus forcing default
> initialization to zero, or should clear() be changed so that it calls
> value_type(0)? In either case, since I am writing a third library which is
> intended to be used with both, my next course of action would be to submit a
> patch to have the behaviour changed.
>
> My vote would be to change uBLAS so that it tests a zero_initialization trait
> at compile time before calling value_type() [implicit initialization to zero],
> otherwise tests a zero_constructible trait before calling value_type(0)
> [explicit initialization to zero], and only otherwise calls value_type()
> [default initialization]. In other words, clear() should try to clear to zero,
> implicitly if possible and explictly where necessary and possible. My GluCat
> library could then define the missing traits for dd_real and qd_real to force
> uBLAS to do the right thing.
>
> The alternative would be to submit a patch to qd to ask that dd_real and
> qd_real be changed to do default-initialization to zero.
> Best, Paul
>
> PS:
>
> I created the file clear_test.cpp with the following contents,
>
> #include <qd/qd_real.h>
>
> #include <boost/numeric/ublas/fwd.hpp>
> #include <boost/numeric/ublas/matrix.hpp>
> #include <boost/numeric/ublas/io.hpp>
>
> template < typename scalar_t >
> void
> test_clear(int dim)
> {
> namespace ublas = boost::numeric::ublas;
> typedef ublas::matrix< scalar_t, ublas::row_major > lhs_t;
>
> std::cout << "default constructor == " << scalar_t() << std::endl;
>
> lhs_t lhs(dim, dim);
> lhs.clear();
>
> std::cout << "lhs:" << lhs << std::endl;
> }
>
> int main(int argc, char ** argv)
> {
> std::cout << "float:" << std::endl;
> test_clear< float >(1);
> std::cout << std::endl;
>
> std::cout << "double:" << std::endl;
> test_clear< double >(1);
> std::cout << std::endl;
>
> unsigned int old_fpu_control;
> fpu_fix_start(&old_fpu_control);
>
> std::cout << "dd_real:" << std::endl;
> test_clear< dd_real >(1);
> std::cout << std::endl;
>
> std::cout << "qd_real:" << std::endl;
> test_clear< qd_real >(1);
> std::cout << std::endl;
>
> return 0;
> }
>
> It produces the following output:
>
> ./clear_test.exe
> float:
> default constructor == 0
> lhs:[1,1]((0))
>
> double:
> default constructor == 0
> lhs:[1,1]((0))
>
> dd_real:
> default constructor == ERROR (dd_real::to_digits): can't compute exponent.
> .e+00
> lhs:ERROR (dd_real::to_digits): can't compute exponent.
> [1,1]((h. <strange characters>
>
> qd_real:
> default constructor == 2.083705e-309
> terminate called after throwing an instance of 'std::logic_error'
> what(): (qd_real::to_digits): can't compute exponent.
> lhs:Aborted
>
> After changing all instances of clear() in matrix.h to use value_type(0)
> rather than value_type(), and changing the "default constructor" output line
> in clear_test.cpp to:
>
> std::cout << "zero constructor == " << scalar_t(0) << std::endl;
>
> the output becomes:
>
> ./clear_test.exe
> float:
> zero constructor == 0
> lhs:[1,1]((0))
>
> double:
> zero constructor == 0
> lhs:[1,1]((0))
>
> dd_real:
> zero constructor == 0.000000e+00
> lhs:[1,1]((0.000000e+00))
>
> qd_real:
> zero constructor == 0.000000e+00
> lhs:[1,1]((0.000000e+00))
>
> _______________________________________________
> ublas mailing list
> ublas_at_[hidden]
> http://lists.boost.org/mailman/listinfo.cgi/ublas
> Sent to: nasos_i_at_[hidden]
_________________________________________________________________
Hotmail has tools for the New Busy. Search, chat and e-mail from your inbox.
http://www.windowslive.com/campaign/thenewbusy?ocid=PID27925::T:WLMTAGL:ON:WL:en-US:WM_HMP:032010_1