[Trilinos-Users] Trilinos 10.0: Tpetra issues
Heroux, Michael A
maherou at sandia.gov
Wed Oct 7 16:29:35 MDT 2009
Nico,
The first issue you experienced is a bug. I have fixed it and it will go out with the 10.0.1 release. I have included a patch file with the changes.
The second issue is a change in semantics from Epetra and Tpetra. Tpetra supports device compute capabilities where memory is not necessarily directly accessible by reference. As a result, you must add one line of code:
Teuchos::ArrayRCP<const double> w_view = w.get1dView();
and then use w_view instead of w in your code.
This process allows us to safely and efficiently "checkout" the data associated with a vector (or multivector). You can then use operator[] with w_view. When the w_view variable goes out of scope, or you reassign its value, the view of w is released. There are also non-const versions of these methods that allow you to change the values (all of these methods are part of Tpetra::MultiVector). If the data is accessible by reference it will already be updated. For devices, when the view is released, the original values will be updated automatically.
I have included a new version of your code that, with the patch applied, works for me.
Please let me know if there are any more issues.
Mike
On 10/6/09 3:42 PM, "Schlömer Nico" <nico.schloemer at ua.ac.be> wrote:
Hi,
two (short?) questions on Tpetra on Trilinos 10.0 here:
I'm just modifying the code from using 10.0 rather than the June alpha release, and it seems that for Tpetra::MultiVector's the handy [] operators have disappeared; same thing for Tpetra::Vector's.
How can I access elements now? I understand that getVector(k) returns the kth vector, but then I'm stuck.
What I try to do is merely print the contents of one vector in its full glory to std::cout, would there possibly be a smarter way to do so than to loop over the the vector's components?
Second, I can't seem to get going with Tpetra::Vector easily. Something that fails to compile is
Tpetra::Vector<double,int> w = Tpetra::Vector<double,int>( myMap );
w.replaceGlobalValue( 0, 0.333 );
with
======================= *snip* =======================
/usr/local/trilinos/10.0/include/Tpetra_Vector.hpp: In member function 'void Tpetra::Vector<S, LO, GO, N>::replaceGlobalValue(GlobalOrdinal, const Scalar&) [with Scalar = double, LocalOrdinal = int, GlobalOrdinal = int, Node = Kokkos::SerialNode]':
tpetraTest.cpp:25: instantiated from here
/usr/local/trilinos/10.0/include/Tpetra_Vector.hpp:67: error: no matching function for call to 'Tpetra::Vector<double, int, int, Kokkos::SerialNode>::replaceGlobalValue(int&, int, const double&)'
/usr/local/trilinos/10.0/include/Tpetra_Vector.hpp:66: note: candidates are: void Tpetra::Vector<S, LO, GO, N>::replaceGlobalValue(GlobalOrdinal, const Scalar&) [with Scalar = double, LocalOrdinal = int, GlobalOrdinal = int, Node = Kokkos::SerialNode]
make: *** [tpetraTest.o] Error 1
======================= *snap* =======================
(where myMap is of course previously defined), and I totally don't understand why. Replacing the Tpetra::Vector's by regular Tpetra::MultiVector's (with one component) can be used as a work-around for that, though.
Any help appreciated!
Nico
-------------- next part --------------
An HTML attachment was scrubbed...
URL: https://software.sandia.gov/pipermail/trilinos-users/attachments/20091007/ce3717d4/attachment.html
-------------- next part --------------
A non-text attachment was scrubbed...
Name: tpetraTest.cpp
Type: application/octet-stream
Size: 1330 bytes
Desc: tpetraTest.cpp
Url : https://software.sandia.gov/pipermail/trilinos-users/attachments/20091007/ce3717d4/attachment.obj
-------------- next part --------------
A non-text attachment was scrubbed...
Name: patchfile
Type: application/octet-stream
Size: 640 bytes
Desc: patchfile
Url : https://software.sandia.gov/pipermail/trilinos-users/attachments/20091007/ce3717d4/attachment-0001.obj
More information about the Trilinos-Users
mailing list