Support Software for Vector Reduction/Transformation Operators  Version of the Day
RTOpPack_SPMD_apply_op_def.hpp
1 // @HEADER
2 // ***********************************************************************
3 //
4 // RTOp: Interfaces and Support Software for Vector Reduction Transformation
5 // Operations
6 // Copyright (2006) Sandia Corporation
7 //
8 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
9 // license for use of this work by or on behalf of the U.S. Government.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov)
39 //
40 // ***********************************************************************
41 // @HEADER
42 
43 #ifndef RTOPPACK_SPMD_APPLY_OP_DEF_HPP
44 #define RTOPPACK_SPMD_APPLY_OP_DEF_HPP
45 
46 #include "RTOpPack_SPMD_apply_op_decl.hpp"
47 #include "Teuchos_Workspace.hpp"
48 #include "Teuchos_CommHelpers.hpp"
49 
50 
51 //
52 // Implementation-only utlities!
53 //
54 
55 
56 namespace RTOpPack {
57 
58 
59 RCP<FancyOStream>& spmdApplyOpDumpOut();
60 
61 
62 template<class Scalar>
63 void print( const ConstSubVectorView<Scalar> &v, Teuchos::FancyOStream &out_arg )
64 {
66  Teuchos::OSTab tab(out);
67  *out << "globalOffset="<<v.globalOffset()<<"\n";
68  *out << "subDim="<<v.subDim()<<"\n";
69  *out << "values:\n";
70  tab.incrTab();
71  for( int i = 0; i < v.subDim(); ++i )
72  *out << " " << v(i) << ":" << (v.globalOffset()+i);
73  *out << "\n";
74 }
75 
76 
77 } // namespace RTOpPack
78 
79 
80 // ///////////////////////////
81 // Template implementations
82 
83 
84 //
85 // Misc Helper functions
86 //
87 
88 
89 template<class PrimitiveScalar>
90 int RTOpPack::serializedSize(
91  int num_values,
92  int num_indexes,
93  int num_chars
94  )
95 {
96  return 3 * sizeof(index_type)
97  + num_values * sizeof(PrimitiveScalar)
98  + num_indexes * sizeof(index_type)
99  + num_chars * sizeof(char_type);
100 }
101 
102 
103 template<class Scalar>
104 void RTOpPack::serialize(
105  const RTOpT<Scalar> &op,
106  Ordinal num_values,
107  Ordinal num_indexes,
108  Ordinal num_chars,
109  const ReductTarget &reduct_obj,
110  char reduct_obj_ext[]
111  )
112 {
113  using Teuchos::arrayView;
114  typedef typename RTOpT<Scalar>::primitive_value_type primitive_value_type;
118  const Ordinal
119  prim_value_type_size = PVTST::fromCountToDirectBytes(1),
120  index_type_size = ITST::fromCountToDirectBytes(1);
121  //char_type_size = CTST::fromCountToDirectBytes(1);
122  const Ordinal
123  num_values_off = 0,
124  num_indexes_off = num_values_off + index_type_size,
125  num_chars_off = num_indexes_off + index_type_size,
126  values_off = num_chars_off + index_type_size,
127  indexes_off = values_off + num_values * prim_value_type_size,
128  chars_off = indexes_off + num_indexes * index_type_size;
129  ITST::serialize(1, &num_values, index_type_size, &reduct_obj_ext[num_values_off]);
130  ITST::serialize(1, &num_indexes, index_type_size, &reduct_obj_ext[num_indexes_off]);
131  ITST::serialize(1, &num_chars, index_type_size, &reduct_obj_ext[num_chars_off]);
132  op.extract_reduct_obj_state(
133  reduct_obj,
134  arrayView(PVTST::convertFromCharPtr(&reduct_obj_ext[values_off]), num_values),
135  arrayView(ITST::convertFromCharPtr(&reduct_obj_ext[indexes_off]), num_indexes),
136  arrayView(CTST::convertFromCharPtr(&reduct_obj_ext[chars_off]), num_chars)
137  );
138  // ToDo: Change above implementation to only require indirect serialization!
139 }
140 
141 
142 template<class Scalar>
143 void RTOpPack::deserialize(
144  const RTOpT<Scalar> &op,
145  int num_values_in,
146  int num_indexes_in,
147  int num_chars_in,
148  const char reduct_obj_ext[],
149  ReductTarget *reduct_obj
150  )
151 {
152  using Teuchos::arrayView;
153  typedef typename RTOpT<Scalar>::primitive_value_type primitive_value_type;
157  const Ordinal
158  prim_value_type_size = PVTST::fromCountToDirectBytes(1),
159  index_type_size = ITST::fromCountToDirectBytes(1);
160  //char_type_size = CTST::fromCountToDirectBytes(1);
161  const Ordinal
162  num_values_off = 0,
163  num_indexes_off = num_values_off + index_type_size,
164  num_chars_off = num_indexes_off + index_type_size,
165  values_off = num_chars_off + index_type_size,
166  indexes_off = values_off + num_values_in * prim_value_type_size,
167  chars_off = indexes_off + num_indexes_in * index_type_size;
168 #ifdef RTOP_DEBUG
169  Ordinal num_values = -1, num_indexes = -1, num_chars = -1;
170  ITST::deserialize(index_type_size, &reduct_obj_ext[num_values_off], 1, &num_values);
171  ITST::deserialize(index_type_size, &reduct_obj_ext[num_indexes_off], 1, &num_indexes);
172  ITST::deserialize(index_type_size, &reduct_obj_ext[num_chars_off], 1, &num_chars);
174  !(
175  num_values==num_values_in && num_indexes==num_indexes_in
176  && num_chars==num_chars_in ),
177  std::logic_error,
178  "Error: RTOp="<<op.op_name()
179  << ", num_values="<<num_values<<", num_values_in="<<num_values_in
180  << ", num_indexes="<<num_indexes<<", num_indexes_in="<<num_indexes_in
181  << ", num_chars="<<num_chars<<", num_chars_in="<<num_chars_in
182  );
183 #endif
184  op.load_reduct_obj_state(
185  arrayView(PVTST::convertFromCharPtr(&reduct_obj_ext[values_off]), num_values_in),
186  arrayView(ITST::convertFromCharPtr(&reduct_obj_ext[indexes_off]), num_indexes_in),
187  arrayView(CTST::convertFromCharPtr(&reduct_obj_ext[chars_off]), num_chars_in),
188  Teuchos::ptr(reduct_obj)
189  );
190  // ToDo: Change above implementation to only require indirect serialization!
191 }
192 
193 
194 namespace RTOpPack {
195 
196 
197 //
198 // ReductTargetSerializer
199 //
200 
201 
202 template<class Scalar>
204  const Teuchos::RCP<const RTOpT<Scalar> > &op
205  )
206  :op_(op.assert_not_null())
207 {
208  using Teuchos::outArg;
209  typedef typename RTOpT<Scalar>::primitive_value_type PrimitiveScalar;
210  op_->get_reduct_type_num_entries(
211  outArg(num_values_), outArg(num_indexes_), outArg(num_chars_) );
212  reduct_obj_ext_size_ =
213  serializedSize<PrimitiveScalar>(num_values_,num_indexes_,num_chars_);
214 }
215 
216 
217 template<class Scalar>
218 index_type
220 {
221  return reduct_obj_ext_size_ * count;
222 }
223 
224 
225 template<class Scalar>
227  const index_type count
228  ,const ReductTarget * const reduct_objs[]
229  ,const index_type bytes
230  ,char charBuffer[]
231  ) const
232 {
233 #ifdef RTOP_DEBUG
234  TEUCHOS_TEST_FOR_EXCEPT( !(count > 0) );
235  TEUCHOS_TEST_FOR_EXCEPT( !reduct_objs );
236  TEUCHOS_TEST_FOR_EXCEPT( !(bytes==this->getBufferSize(count)) );
237  TEUCHOS_TEST_FOR_EXCEPT( !charBuffer );
238 #endif
239  Ordinal offset = 0;
240  for( Ordinal i = 0; i < count; ++i, offset += reduct_obj_ext_size_ ) {
241  RTOpPack::serialize(
242  *op_,num_values_,num_indexes_,num_chars_
243  ,*reduct_objs[i],&charBuffer[offset]
244  );
245  }
246 }
247 
248 
249 template<class Scalar>
252 {
253  return op_->reduct_obj_create();
254 }
255 
256 template<class Scalar>
258  const index_type bytes
259  ,const char charBuffer[]
260  ,const index_type count
261  ,ReductTarget * const reduct_objs[]
262  ) const
263 {
264 #ifdef RTOP_DEBUG
265  TEUCHOS_TEST_FOR_EXCEPT( !(bytes > 0) );
266  TEUCHOS_TEST_FOR_EXCEPT( !charBuffer );
267  TEUCHOS_TEST_FOR_EXCEPT( !(bytes==getBufferSize(count)) );
268  TEUCHOS_TEST_FOR_EXCEPT( !reduct_objs );
269 #endif
270  Ordinal offset = 0;
271  for( Ordinal i = 0; i < count; ++i, offset += reduct_obj_ext_size_ ) {
272  RTOpPack::deserialize(
273  *op_,num_values_,num_indexes_,num_chars_
274  ,&charBuffer[offset],reduct_objs[i]
275  );
276  }
277 }
278 
279 
280 //
281 // ReductTargetReductionOp
282 //
283 
284 
285 template<class Scalar>
287  const Teuchos::RCP<const RTOpT<Scalar> > &op
288  )
289  :op_(op)
290 {}
291 
292 
293 template<class Scalar>
295  const Ordinal count,
296  const ReductTarget*const inBuffer[],
297  ReductTarget*const inoutBuffer[]
298  ) const
299 {
300  for( Ordinal i = 0; i < count; ++i )
301  op_->reduce_reduct_objs( *inBuffer[i], Teuchos::ptr(inoutBuffer[i]) );
302 }
303 
304 
305 } // namespace RTOpPack
306 
307 
308 template<class Scalar>
309 void RTOpPack::SPMD_all_reduce(
310  const Teuchos::Comm<index_type> *comm,
311  const RTOpT<Scalar> &op,
312  const int num_cols,
313  const ReductTarget*const i_reduct_objs[],
314  ReductTarget*const reduct_objs[]
315  )
316 {
317  using Teuchos::Workspace;
318  using Teuchos::reduceAll;
320  Workspace<Teuchos::RCP<ReductTarget> >
321  i_i_reduct_objs( wss, num_cols );
322  Workspace<ReductTarget*>
323  _i_i_reduct_objs( wss, num_cols );
324  for( int kc = 0; kc < num_cols; ++kc ) {
325  i_i_reduct_objs[kc] = op.reduct_obj_create();
326  _i_i_reduct_objs[kc] = &*i_i_reduct_objs[kc];
327  }
328  ReductTargetSerializer<Scalar>
329  serializer(Teuchos::rcpFromRef(op));
330  ReductTargetReductionOp<Scalar>
331  reductOp(Teuchos::rcpFromRef(op));
332  reduceAll<Ordinal>(
333  *comm, serializer, reductOp,
334  num_cols, &i_reduct_objs[0], &_i_i_reduct_objs[0]);
335  for( int kc = 0; kc < num_cols; ++kc ) {
336  op.reduce_reduct_objs(*_i_i_reduct_objs[kc], Teuchos::ptr(reduct_objs[kc]));
337  }
338 }
339 
340 
341 template<class Scalar>
342 void RTOpPack::SPMD_apply_op(
343  const Teuchos::Comm<index_type> *comm,
344  const RTOpT<Scalar> &op,
345  const int num_vecs,
346  const RTOpPack::ConstSubVectorView<Scalar> sub_vecs[],
347  const int num_targ_vecs,
348  const RTOpPack::SubVectorView<Scalar> targ_sub_vecs[],
349  ReductTarget *reduct_obj
350  )
351 {
352  ReductTarget* reduct_objs[] = { reduct_obj };
353  SPMD_apply_op(
354  comm,op,1,num_vecs,sub_vecs,num_targ_vecs,targ_sub_vecs
355  ,reduct_obj ? reduct_objs : NULL
356  );
357 }
358 
359 
361 template<class Scalar>
362 void RTOpPack::SPMD_apply_op(
363  const Teuchos::Comm<index_type> *comm,
364  const RTOpT<Scalar> &op,
365  const int num_cols,
366  const int num_multi_vecs,
367  const RTOpPack::ConstSubMultiVectorView<Scalar> sub_multi_vecs[],
368  const int num_targ_multi_vecs,
369  const RTOpPack::SubMultiVectorView<Scalar> targ_sub_multi_vecs[],
370  RTOpPack::ReductTarget*const reduct_objs[]
371  )
372 {
373  using Teuchos::arcp;
374  using Teuchos::Workspace;
376  int k, j, off;
377  Workspace<ConstSubVectorView<Scalar> > c_sub_vecs(wss,num_multi_vecs*num_cols);
378  if(sub_multi_vecs) {
379  for( off = 0, j = 0; j < num_cols; ++j ) {
380  for( k = 0; k < num_multi_vecs; ++k ) {
381  const ConstSubMultiVectorView<Scalar> &mv = sub_multi_vecs[k];
382  if (mv.subDim()) {
383  c_sub_vecs[off++].initialize(mv.globalOffset(), mv.subDim(),
384  arcp(&mv(0,j), 0, mv.subDim(), false), 1);
385  }
386  else {
387  c_sub_vecs[off++].initialize(mv.globalOffset(), mv.subDim(),
388  Teuchos::null, 1);
389  }
390  }
391  }
392  }
393  Workspace<SubVectorView<Scalar> > c_targ_sub_vecs(wss,num_targ_multi_vecs*num_cols);
394  if(targ_sub_multi_vecs) {
395  for( off = 0, j = 0; j < num_cols; ++j ) {
396  for( k = 0; k < num_targ_multi_vecs; ++k ) {
397  const SubMultiVectorView<Scalar> &mv = targ_sub_multi_vecs[k];
398  ArrayRCP<Scalar> mv_j = Teuchos::null;
399  if (mv.subDim()) { mv_j = arcp(&mv(0,j), 0, mv.subDim(), false); }
400  c_targ_sub_vecs[off++].initialize(mv.globalOffset(), mv.subDim(), mv_j, 1);
401  }
402  }
403  }
404  SPMD_apply_op(
405  comm,op,num_cols
406  ,num_multi_vecs, num_multi_vecs && sub_multi_vecs ? &c_sub_vecs[0] : NULL
407  ,num_targ_multi_vecs, num_targ_multi_vecs && targ_sub_multi_vecs ? &c_targ_sub_vecs[0] : NULL
408  ,reduct_objs
409  );
410 }
411 
412 
413 template<class Scalar>
414 void RTOpPack::SPMD_apply_op(
415  const Teuchos::Comm<index_type> *comm,
416  const RTOpT<Scalar> &op,
417  const int num_cols,
418  const int num_vecs,
419  const ConstSubVectorView<Scalar> sub_vecs[],
420  const int num_targ_vecs,
421  const SubVectorView<Scalar> sub_targ_vecs[],
422  ReductTarget*const reduct_objs[]
423  )
424 {
425  using Teuchos::arrayView;
426  Teuchos::RCP<Teuchos::FancyOStream> out = spmdApplyOpDumpOut();
427  Teuchos::OSTab tab(out);
428  if (nonnull(out)) {
429  *out << "\nEntering RTOpPack::SPMD_apply_op(...) ...\n";
430  *out
431  << "\ncomm = " << (comm?comm->description():"NULL")
432  << "\nop = " << op.description()
433  << "\nnum_cols = " << num_cols
434  << "\nnum_vecs = " << num_vecs
435  << "\nnum_targ_vecs = " << num_targ_vecs
436  << "\n";
437  if( num_vecs && sub_vecs ) {
438  *out << "\nInput vectors:\n";
439  Teuchos::OSTab tab2(out);
440  for( int kc = 0; kc < num_cols; ++kc ) {
441  for( int k = 0; k < num_vecs; ++k ) {
442  *out << "\nvecs["<<kc<<","<<k<<"] =\n";
443  print(sub_vecs[kc*num_vecs+k],*out);
444  }
445  }
446  }
447  if( num_targ_vecs && sub_targ_vecs ) {
448  *out << "\nInput/output vectors *before* transforamtion:\n";
449  Teuchos::OSTab tab2(out);
450  for( int kc = 0; kc < num_cols; ++kc ) {
451  for( int k = 0; k < num_targ_vecs; ++k ) {
452  *out << "\nvecs["<<kc<<","<<k<<"] =\n";
453  print(sub_targ_vecs[kc*num_targ_vecs+k],*out);
454  }
455  }
456  }
457  if(reduct_objs) {
458  *out << "\nInput/output reduction objects *before* reduction:\n";
459  Teuchos::OSTab tab2(out);
460  for( int kc = 0; kc < num_cols; ++kc ) {
461  *out
462  << "\nreduct_objs["<<kc<<"] =\n"
463  << describe(*reduct_objs[kc],Teuchos::VERB_EXTREME);
464  }
465  }
466  }
467  using Teuchos::Workspace;
469  if( reduct_objs == NULL && sub_vecs == NULL && sub_targ_vecs == NULL ) {
470  // This is a transformation operation with no data on this processor.
471  // Therefore, we can just exist!
472  }
473  else {
474  const int localSubDim =
475  ( num_vecs
476  ? ( sub_vecs ? sub_vecs[0].subDim() : 0 )
477  : ( sub_targ_vecs ? sub_targ_vecs[0].subDim() : 0 )
478  );
479  // See if we need to do any global communication at all?
480  if( comm==NULL || reduct_objs == NULL ) {
481  if( ( sub_vecs || sub_targ_vecs ) && localSubDim ) {
482  for( int kc = 0; kc < num_cols; ++kc ) {
483  op.apply_op(
484  arrayView(sub_vecs+kc*num_vecs, num_vecs),
485  arrayView(sub_targ_vecs+kc*num_targ_vecs, num_targ_vecs),
486  reduct_objs ? Teuchos::ptr(reduct_objs[kc]) : Teuchos::null
487  );
488  }
489  }
490  }
491  else {
492  // Check the preconditions for excluding empty target vectors.
494  ( ( num_vecs && !sub_vecs) || ( num_targ_vecs && !sub_targ_vecs) ) && !( !sub_vecs && !sub_targ_vecs )
495  ,std::logic_error
496  ,"SPMD_apply_op(...): Error, invalid arguments num_vecs = " << num_vecs
497  << ", sub_vecs = " << sub_vecs << ", num_targ_vecs = " << num_targ_vecs
498  << ", sub_targ_vecs = " << sub_targ_vecs
499  );
500  //
501  // There is a non-null reduction target object and we are using
502  // SPMD so we need to reduce it across processors
503  //
504  // Allocate the intermediate target object and perform the
505  // reduction for the vector elements on this processor.
506  //
507  Workspace<Teuchos::RCP<ReductTarget> >
508  i_reduct_objs( wss, num_cols );
509  for( int kc = 0; kc < num_cols; ++kc ) {
510  i_reduct_objs[kc] = op.reduct_obj_create();
511  if( ( sub_vecs || sub_targ_vecs ) && localSubDim ) {
512  op.apply_op(
513  arrayView(sub_vecs+kc*num_vecs, num_vecs),
514  arrayView(sub_targ_vecs+kc*num_targ_vecs, num_targ_vecs),
515  i_reduct_objs[kc].ptr()
516  );
517  }
518  }
519  if(nonnull(out)) {
520  if(reduct_objs) {
521  *out << "\nIntermediate reduction objects in this process before global reduction:\n";
522  Teuchos::OSTab tab2(out);
523  for( int kc = 0; kc < num_cols; ++kc ) {
524  *out
525  << "\ni_reduct_objs["<<kc<<"] =\n"
526  << describe(*i_reduct_objs[kc],Teuchos::VERB_EXTREME);
527  }
528  }
529  }
530  //
531  // Reduce the local intermediate reduction objects into the global reduction objects
532  //
533  Workspace<const ReductTarget*>
534  _i_reduct_objs( wss, num_cols );
535  for( int kc = 0; kc < num_cols; ++kc ) {
536  _i_reduct_objs[kc] = &*i_reduct_objs[kc];
537  }
538  if(nonnull(out)) {
539  if(reduct_objs) {
540  *out << "\nPerforming global reduction ...\n";
541  }
542  }
543  SPMD_all_reduce(comm,op,num_cols,&_i_reduct_objs[0],reduct_objs);
544  }
545  }
546  if(nonnull(out)) {
547  if( num_targ_vecs && sub_targ_vecs ) {
548  *out << "\nInput/output vectors *after* transforamtion:\n";
549  Teuchos::OSTab tab2(out);
550  for( int kc = 0; kc < num_cols; ++kc ) {
551  for( int k = 0; k < num_targ_vecs; ++k ) {
552  *out << "\nvecs["<<kc<<","<<k<<"] =\n";
553  print(sub_targ_vecs[kc*num_targ_vecs+k],*out);
554  }
555  }
556  }
557  if(reduct_objs) {
558  *out << "\nInput/output reduction objects *after* reduction:\n";
559  Teuchos::OSTab tab2(out);
560  for( int kc = 0; kc < num_cols; ++kc ) {
561  *out
562  << "\nreduct_objs["<<kc<<"] =\n"
563  << describe(*reduct_objs[kc],Teuchos::VERB_EXTREME);
564  }
565  }
566  *out << "\nLeaving RTOpPack::SPMD_apply_op(...) ...\n";
567  *out << std::flush;
568  }
569 }
570 
571 
572 //
573 // Explicit Template Instaniation Macros
574 //
575 
576 
577 #define RTOPPACK_SPMD_APPLY_OP_INSTANT_SCALAR(SCALAR) \
578  \
579  template int serializedSize<SCALAR >( \
580  int num_values, \
581  int num_indexes, \
582  int num_chars \
583  ); \
584  \
585  template void serialize<SCALAR >( \
586  const RTOpT<SCALAR > &op, \
587  Ordinal num_values, \
588  Ordinal num_indexes, \
589  Ordinal num_chars, \
590  const ReductTarget &reduct_obj, \
591  char reduct_obj_ext[] \
592  ); \
593  \
594  template void deserialize<SCALAR >( \
595  const RTOpT<SCALAR > &op, \
596  int num_values_in, \
597  int num_indexes_in, \
598  int num_chars_in, \
599  const char reduct_obj_ext[], \
600  ReductTarget *reduct_obj \
601  ); \
602  \
603  template class ReductTargetSerializer<SCALAR >; \
604  \
605  template class ReductTargetReductionOp<SCALAR >; \
606  \
607  template void SPMD_all_reduce<SCALAR >( \
608  const Teuchos::Comm<index_type> *comm, \
609  const RTOpT<SCALAR > &op, \
610  const int num_cols, \
611  const ReductTarget*const i_reduct_objs[], \
612  ReductTarget*const reduct_objs[] \
613  ); \
614  \
615  template void SPMD_apply_op<SCALAR >( \
616  const Teuchos::Comm<index_type> *comm, \
617  const RTOpT<SCALAR > &op, \
618  const int num_vecs, \
619  const RTOpPack::ConstSubVectorView<SCALAR > sub_vecs[], \
620  const int num_targ_vecs, \
621  const RTOpPack::SubVectorView<SCALAR > targ_sub_vecs[], \
622  ReductTarget *reduct_obj \
623  ); \
624  \
625  template void SPMD_apply_op<SCALAR >( \
626  const Teuchos::Comm<index_type> *comm, \
627  const RTOpT<SCALAR > &op, \
628  const int num_cols, \
629  const int num_multi_vecs, \
630  const RTOpPack::ConstSubMultiVectorView<SCALAR > sub_multi_vecs[], \
631  const int num_targ_multi_vecs, \
632  const RTOpPack::SubMultiVectorView<SCALAR > targ_sub_multi_vecs[], \
633  RTOpPack::ReductTarget*const reduct_objs[] \
634  ); \
635  \
636  template void SPMD_apply_op<SCALAR >( \
637  const Teuchos::Comm<index_type> *comm, \
638  const RTOpT<SCALAR > &op, \
639  const int num_cols, \
640  const int num_vecs, \
641  const ConstSubVectorView<SCALAR > sub_vecs[], \
642  const int num_targ_vecs, \
643  const SubVectorView<SCALAR > sub_targ_vecs[], \
644  ReductTarget*const reduct_objs[] \
645  );
646 
647 
648 #endif // RTOPPACK_SPMD_APPLY_OP_DEF_HPP
void deserialize(const index_type bytes, const char charBuffer[], const index_type count, ReductTarget *const reduct_objs[]) const
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
VERB_EXTREME
ReductionOp subclass for ReductTarget objects.
virtual std::string description() const
Serializer subclass for ReductTarget objects.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< ReductTarget > createObj() const
void serialize(const index_type count, const ReductTarget *const reduct_objs[], const index_type bytes, char charBuffer[]) const
void initialize(Ordinal globalOffset_in, Ordinal subDim_in, Ordinal colOffset_in, Ordinal numSubCols_in, const ArrayRCP< Scalar > &values_in, Ordinal leadingDim_in)
PrimitiveTypeTraits< Scalar, Scalar >::primitiveType primitive_value_type
bool nonnull(const boost::shared_ptr< T > &p)
void reduce(const Ordinal count, const ReductTarget *const inBuffer[], ReductTarget *const inoutBuffer[]) const
void initialize(Ordinal globalOffset_in, Ordinal subDim_in, Ordinal colOffset_in, Ordinal numSubCols_in, const ArrayRCP< const Scalar > &values_in, Ordinal leadingDim_in)
index_type getBufferSize(const index_type count) const
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)
TEUCHOSCORE_LIB_DLL_EXPORT Teuchos::RCP< WorkspaceStore > get_default_workspace_store()