43 #ifndef IFPACK2_DETAILS_MULTIVECTORLOCALGATHERSCATTER_HPP
44 #define IFPACK2_DETAILS_MULTIVECTORLOCALGATHERSCATTER_HPP
50 #include "Tpetra_MultiVector.hpp"
51 #include "Tpetra_Map.hpp"
84 template<
class MV_in,
class MV_out>
87 typedef typename MV_in::scalar_type InScalar;
88 typedef typename MV_out::scalar_type OutScalar;
89 typedef typename MV_in::local_ordinal_type LO;
90 typedef typename MV_in::global_ordinal_type GO;
91 typedef typename MV_in::node_type NO;
97 gather (MV_out& X_out,
102 const size_t numRows = X_out.getLocalLength ();
103 const size_t numVecs = X_in.getNumVectors ();
105 for (
size_t j = 0; j < numVecs; ++j) {
106 ArrayRCP<const InScalar> X_in_j = X_in.getData(j);
107 ArrayRCP<OutScalar> X_out_j = X_out.getDataNonConst(j);
108 for (
size_t i = 0; i < numRows; ++i) {
109 const size_t i_perm = perm[i];
110 X_out_j[i] = X_in_j[i_perm];
127 const size_t numBlocks = X_out.getLocalLength() / blockSize;
128 const size_t numVecs = X_in.getNumVectors ();
130 for (
size_t j = 0; j < numVecs; ++j) {
131 ArrayRCP<const InScalar> X_in_j = X_in.getData(j);
132 ArrayRCP<OutScalar> X_out_j = X_out.getDataNonConst(j);
133 for (
size_t i = 0; i < numBlocks; ++i) {
134 const size_t i_perm = perm[i];
135 for (LO k = 0; k < blockSize; k++) {
136 X_out_j[i * blockSize + k] = X_in_j[i_perm * blockSize + k];
144 scatter (MV_in& X_in,
149 const size_t numRows = X_out.getLocalLength();
150 const size_t numVecs = X_in.getNumVectors();
152 for (
size_t j = 0; j < numVecs; ++j) {
153 ArrayRCP<InScalar> X_in_j = X_in.getDataNonConst(j);
154 ArrayRCP<const OutScalar> X_out_j = X_out.getData(j);
155 for (
size_t i = 0; i < numRows; ++i) {
156 const size_t i_perm = perm[i];
157 X_in_j[i_perm] = X_out_j[i];
171 const size_t numBlocks = X_out.getLocalLength() / blockSize;
172 const size_t numVecs = X_in.getNumVectors ();
174 for (
size_t j = 0; j < numVecs; ++j) {
175 ArrayRCP<const InScalar> X_in_j = X_in.getData(j);
176 ArrayRCP<OutScalar> X_out_j = X_out.getDataNonConst(j);
177 for (
size_t i = 0; i < numBlocks; ++i) {
178 const size_t i_perm = perm[i];
179 for (LO k = 0; k < blockSize; k++) {
180 X_in_j[i_perm * blockSize + k] = X_out_j[i * blockSize + k];
190 template<
typename InView,
typename OutView>
191 void gatherViewToView(OutView X_out,
197 for(
size_t j = 0; j < X_out.extent(1); ++j) {
198 for(
size_t i = 0; i < X_out.extent(0); ++i) {
199 const LO i_perm = perm[i];
200 X_out(i, j) = X_in(i_perm, j);
205 template<
typename InView,
typename OutView>
206 void scatterViewToView(InView X_in,
211 for(
size_t j = 0; j < X_out.extent(1); ++j) {
212 for(
size_t i = 0; i < X_out.extent(0); ++i) {
213 const LO i_perm = perm[i];
214 X_in(i_perm, j) = X_out(i, j);
219 template<
typename InView,
typename OutView>
220 void gatherViewToViewBlock(OutView X_out,
227 size_t numBlocks = X_out.extent(0) / blockSize;
228 for(
size_t j = 0; j < X_out.extent(1); ++j) {
229 for(
size_t i = 0; i < numBlocks; ++i) {
230 const LO i_perm = perm[i];
231 for(LO k = 0; k < blockSize; k++) {
232 X_out(i * blockSize + k, j) = X_in(i_perm * blockSize + k, j);
238 template<
typename InView,
typename OutView>
239 void scatterViewToViewBlock(InView X_in,
246 size_t numBlocks = X_out.extent(0) / blockSize;
247 for(
size_t j = 0; j < X_out.extent(1); ++j) {
248 for(
size_t i = 0; i < numBlocks; ++i) {
249 const LO i_perm = perm[i];
250 for(LO k = 0; k < blockSize; k++) {
251 X_in(i_perm * blockSize + k, j) = X_out(i * blockSize + k, j);
260 template<
typename InView>
261 void gatherMVtoView(MV_out X_out,
267 size_t numRows = X_out.getLocalLength();
268 for(
size_t j = 0; j < X_out.getNumVectors(); ++j) {
270 for(
size_t i = 0; i < numRows; ++i) {
271 const LO i_perm = perm[i];
272 X_out_j[i] = X_in(i_perm, j);
277 template<
typename InView>
278 void scatterMVtoView(InView X_in,
282 size_t numRows = X_out.getLocalLength();
284 for(
size_t j = 0; j < X_in.extent(1); ++j) {
286 for(
size_t i = 0; i < numRows; ++i) {
287 const LO i_perm = perm[i];
288 X_in(i_perm, j) = X_out_j[i];
293 template<
typename InView>
294 void gatherMVtoViewBlock(MV_out X_out,
300 size_t numBlocks = X_out.getLocalLength() / blockSize;
302 for(
size_t j = 0; j < X_out.getNumVectors(); ++j) {
304 for(
size_t i = 0; i < numBlocks; ++i) {
305 const LO i_perm = perm[i];
306 for(LO k = 0; k < blockSize; k++) {
307 X_out_j[i * blockSize + k] = X_in(i_perm * blockSize + k, j);
313 template<
typename InView>
314 void scatterMVtoViewBlock(InView X_in,
319 size_t numBlocks = X_out.getLocalLength() / blockSize;
321 for(
size_t j = 0; j < X_in.extent(1); ++j) {
323 for(
size_t i = 0; i < numBlocks; ++i) {
324 const LO i_perm = perm[i];
325 for(LO k = 0; k < blockSize; k++) {
326 X_in(i_perm * blockSize + k, j) = X_out_j[i * blockSize + k];
336 #endif // IFPACK2_DETAILS_MULTIVECTORLOCALGATHERSCATTER_HPP
Implementation detail of Ifpack2::Container subclasses.
Definition: Ifpack2_Details_MultiVectorLocalGatherScatter.hpp:85