Hi Karl,

Karl Rupp <r...@iue.tuwien.ac.at> writes:
> It seems like you found a bug. In the iterative solvers we have the lines:
>        VectorType result = rhs;
>        viennacl::traits::clear(result);
>
> If VectorType is viennacl::vector, everything is okay. However, if 
> VectorType is viennacl::vector_base, then a missing copy-CTOR results in 
> a compiler-generated, 'flat' copy. The subsequent clear() thus operates 
> on the memory handle copied from 'rhs', clearing the data.
> Ticket created: https://github.com/viennacl/viennacl-dev/issues/60

That would explain why the vector mysteriously disappears. Notably, when
I put the print commands in my C++ code, I put them *after* the solver
call (probably a mistake to put them there, in hindsight, but
nonetheless it seems to have been useful!).

Why does *_base need a copy constructor anyway? Why not force copies to
be made using derived types? It seems strange to have copies of _base
instances lying around without the higher-level objects to which they
normally refer.

And this explains why the direct solvers work: they have different
prototypes. Instead of VectorType, you have a return type of
vector<NumericT>, which of course skips this bug. Since it's weird to
have a solver returning a vector_base object in any case, why not just
return vector<T> for the iterative solvers, too?

In fact, the different styles of prototype were a bit awkward for me:
g++ was getting confused, forcing me to split the code into three files
(direct_solvers.cpp, iterative_solvers.cpp and solve_op_func.hpp, where
solve_op_func.hpp contained the call that g++ previously thought
"ambiguous"). It turns out that if you include both direct_solve.hpp and
gmres.hpp into your program, g++ 4.8 can't disambiguate the call based
on the tag class, as you'd expect.

>> Does anyone know why that might be? This doesn't seem to happen
>> elsewhere.. (Notably, for instance, the direct solvers work, and those
>> calls use matrix_base and vector_base prototypes).
>
> Hmm, seems like this requires some temporary to fix things up. Could you 
> please create a temporary object rhs_copy of type viennacl_vector<> and 
> pass that to the iterative solvers? For example, replace
>
>   result = viennacl::solve(A, b, cg_tag());
>
> by
>   viennacl::vector<T> b_copy(b);
>   result = viennacl::solve(A, b_copy, cg_tag());
>
> The additional temporary object won't matter for overall performance.

I'll do that. And I'll create the copy in only those cases where the
object isn't exactly of vector<T> type (for instance, where we have a
proxy instead). I expect those cases to be fairly rare in real usage, of
course.


Thanks for your help!

Toby


------------------------------------------------------------------------------
Managing the Performance of Cloud-Based Applications
Take advantage of what the Cloud has to offer - Avoid Common Pitfalls.
Read the Whitepaper.
http://pubads.g.doubleclick.net/gampad/clk?id=121054471&iu=/4140/ostg.clktrk
_______________________________________________
ViennaCL-devel mailing list
ViennaCL-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/viennacl-devel

Reply via email to