Thyra_SpmdVectorBase.hpp

Go to the documentation of this file.
00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 //    Thyra: Interfaces and Support for Abstract Numerical Algorithms
00005 //                 Copyright (2004) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 #ifndef THYRA_SPMD_VECTOR_BASE_HPP
00030 #define THYRA_SPMD_VECTOR_BASE_HPP
00031 
00032 #include "Thyra_SpmdVectorBaseDecl.hpp"
00033 #include "Thyra_VectorDefaultBase.hpp"
00034 #include "Thyra_SpmdVectorSpaceDefaultBase.hpp"
00035 #include "Thyra_apply_op_helper.hpp"
00036 #include "RTOp_parallel_helpers.h"
00037 #include "RTOpPack_SPMD_apply_op.hpp"
00038 #include "Teuchos_Workspace.hpp"
00039 #include "Teuchos_TestForException.hpp"
00040 #include "Teuchos_dyn_cast.hpp"
00041 
00042 #ifdef THYRA_SPMD_VECTOR_BASE_DUMP
00043 #  include "Teuchos_VerboseObject.hpp"
00044 #endif // THYRA_SPMD_VECTOR_BASE_DUMP
00045 
00046 namespace Thyra {
00047 
00048 template<class Scalar>
00049 SpmdVectorBase<Scalar>::SpmdVectorBase()
00050   :in_applyOp_(false)
00051   ,globalDim_(0)
00052   ,localOffset_(-1)
00053   ,localSubDim_(0)
00054 {}
00055 
00056 // Virtual methods with default implementations
00057 
00058 template<class Scalar>
00059 void SpmdVectorBase<Scalar>::getLocalData( const Scalar** values, Index* stride ) const
00060 {
00061   const_cast<SpmdVectorBase<Scalar>*>(this)->getLocalData(const_cast<Scalar**>(values),stride);
00062 }
00063 
00064 template<class Scalar>
00065 void SpmdVectorBase<Scalar>::freeLocalData( const Scalar* values ) const
00066 {
00067   const_cast<SpmdVectorBase<Scalar>*>(this)->commitLocalData(const_cast<Scalar*>(values));
00068 }
00069 
00070 template<class Scalar>
00071 void SpmdVectorBase<Scalar>::applyOp(
00072   const Teuchos::Comm<Index>      *comm_in
00073   ,const RTOpPack::RTOpT<Scalar>  &op
00074   ,const int                      num_vecs
00075   ,const VectorBase<Scalar>*const vecs[]
00076   ,const int                      num_targ_vecs
00077   ,VectorBase<Scalar>*const       targ_vecs[]
00078   ,RTOpPack::ReductTarget         *reduct_obj
00079   ,const Index                    first_ele_offset_in
00080   ,const Index                    sub_dim_in
00081   ,const Index                    global_offset_in
00082   ) const
00083 {
00084 #ifdef THYRA_SPMD_VECTOR_BASE_DUMP
00085   Teuchos::RefCountPtr<Teuchos::FancyOStream>
00086     out = Teuchos::VerboseObjectBase::getDefaultOStream();
00087   Teuchos::OSTab tab(out);
00088   if(show_dump) {
00089     *out << "\nEntering SpmdVectorBase<Scalar>::applyOp(...) ...\n";
00090     *out
00091       << "\nop = " << typeid(op).name()
00092       << "\nnum_vecs = " << num_vecs
00093       << "\nnum_targ_vecs = " << num_targ_vecs
00094       << "\nreduct_obj = " << reduct_obj
00095       << "\nfirst_ele_offset_in = " << first_ele_offset_in
00096       << "\nsub_dim_in = " << sub_dim_in
00097       << "\nglobal_offset_in = " << global_offset_in
00098       << "\n"
00099       ;
00100   }
00101 #endif // THYRA_SPMD_VECTOR_BASE_DUMP
00102   using Teuchos::dyn_cast;
00103   using Teuchos::Workspace;
00104   Teuchos::WorkspaceStore* wss = Teuchos::get_default_workspace_store().get();
00105   const SpmdVectorSpaceBase<Scalar> &spmdSpc = *spmdSpace();
00106 #ifdef TEUCHOS_DEBUG
00107   // ToDo: Validate input!
00108   TEST_FOR_EXCEPTION(
00109     in_applyOp_, std::invalid_argument
00110     ,"SpmdVectorBase<>::applyOp(...): Error, this method is being entered recursively which is a "
00111     "clear sign that one of the methods acquireDetachedView(...), releaseDetachedView(...) or commitDetachedView(...) "
00112     "was not implemented properly!"
00113     );
00114   Thyra::apply_op_validate_input(
00115     "SpmdVectorBase<>::applyOp(...)",*space()
00116     ,op,num_vecs,vecs,num_targ_vecs,targ_vecs,reduct_obj,first_ele_offset_in,sub_dim_in,global_offset_in
00117     );
00118 #endif
00119   Teuchos::RefCountPtr<const Teuchos::Comm<Index> > comm;
00120   if( comm_in != NULL )
00121     comm = Teuchos::rcp(comm_in,false);
00122   else
00123     comm = spmdSpc.getComm();
00124   // Flag that we are in applyOp()
00125   in_applyOp_ = true;
00126   // First see if this is a locally replicated vector in which case
00127   // we treat this as a local operation only.
00128   const bool locallyReplicated = ( comm_in == NULL && localSubDim_ == globalDim_ );
00129   // Get the overlap in the current process with the input logical sub-vector
00130   // from (first_ele_offset_in,sub_dim_in,global_offset_in)
00131   Index  overlap_first_local_ele_off  = 0;
00132   Index  overlap_local_sub_dim        = 0;
00133   Index  overlap_global_off           = 0;
00134   if(localSubDim_) {
00135     RTOp_parallel_calc_overlap(
00136       globalDim_, localSubDim_, localOffset_, first_ele_offset_in, sub_dim_in, global_offset_in
00137       ,&overlap_first_local_ele_off, &overlap_local_sub_dim, &overlap_global_off
00138       );
00139   }
00140   const Range1D local_rng = (
00141     overlap_first_local_ele_off>=0
00142     ? Range1D( localOffset_ + overlap_first_local_ele_off, localOffset_ + overlap_first_local_ele_off + overlap_local_sub_dim - 1 )
00143     : Range1D::Invalid
00144     );
00145 #ifdef THYRA_SPMD_VECTOR_BASE_DUMP
00146   if(show_dump) {
00147     *out
00148       << "\noverlap_first_local_ele_off = " << overlap_first_local_ele_off
00149       << "\noverlap_local_sub_dim = " << overlap_local_sub_dim
00150       << "\noverlap_global_off = " << overlap_global_off
00151       << "\nlocal_rng = ["<<local_rng.lbound()<<","<<local_rng.ubound()<<"]"
00152       << "\n"
00153       ;
00154   }
00155 #endif // THYRA_SPMD_VECTOR_BASE_DUMP
00156   // Create sub-vector views of all of the *participating* local data
00157   Workspace<RTOpPack::ConstSubVectorView<Scalar> > sub_vecs(wss,num_vecs);
00158   Workspace<RTOpPack::SubVectorView<Scalar> > sub_targ_vecs(wss,num_targ_vecs);
00159   if( overlap_first_local_ele_off >= 0 ) {
00160     if(1){for(int k = 0; k < num_vecs; ++k ) {
00161       vecs[k]->acquireDetachedView( local_rng, &sub_vecs[k] );
00162       sub_vecs[k].setGlobalOffset( overlap_global_off );
00163     }}
00164     if(1){for(int k = 0; k < num_targ_vecs; ++k ) {
00165       targ_vecs[k]->acquireDetachedView( local_rng, &sub_targ_vecs[k] );
00166       sub_targ_vecs[k].setGlobalOffset( overlap_global_off );
00167     }}
00168   }
00169   // Apply the RTOp operator object (all processors must participate)
00170   const bool validSubVecs = ( localSubDim_ > 0 && overlap_first_local_ele_off >= 0 );
00171   RTOpPack::SPMD_apply_op(
00172     locallyReplicated ? NULL : comm.get()                       // comm
00173     ,op                                                         // op
00174     ,num_vecs                                                   // num_vecs
00175     ,num_vecs && validSubVecs ? &sub_vecs[0] : NULL             // sub_vecs
00176     ,num_targ_vecs                                              // num_targ_vecs
00177     ,num_targ_vecs && validSubVecs ? &sub_targ_vecs[0] : NULL   // targ_sub_vecs
00178     ,reduct_obj                                                 // reduct_obj
00179     );
00180   // Free and commit the local data
00181   if( overlap_first_local_ele_off >= 0 ) {
00182     if(1){for(int k = 0; k < num_vecs; ++k ) {
00183       sub_vecs[k].setGlobalOffset(local_rng.lbound());
00184       vecs[k]->releaseDetachedView( &sub_vecs[k] );
00185     }}
00186     if(1){for(int k = 0; k < num_targ_vecs; ++k ) {
00187       sub_targ_vecs[k].setGlobalOffset(local_rng.lbound());
00188       targ_vecs[k]->commitDetachedView( &sub_targ_vecs[k] );
00189     }}
00190   }
00191   // Flag that we are leaving applyOp()
00192   in_applyOp_ = false;
00193 #ifdef THYRA_SPMD_VECTOR_BASE_DUMP
00194   if(show_dump) {
00195     *out << "\nLeaving SpmdVectorBase<Scalar>::applyOp(...) ...\n";
00196   }
00197 #endif // THYRA_SPMD_VECTOR_BASE_DUMP
00198 }
00199 
00200 // Overridden from VectorBase
00201 
00202 template<class Scalar>
00203 Teuchos::RefCountPtr<const VectorSpaceBase<Scalar> >
00204 SpmdVectorBase<Scalar>::space() const
00205 {
00206   return spmdSpace();
00207 }
00208 
00209 template<class Scalar>
00210 void SpmdVectorBase<Scalar>::applyOp(
00211   const RTOpPack::RTOpT<Scalar>   &op
00212   ,const int                      num_vecs
00213   ,const VectorBase<Scalar>*const vecs[]
00214   ,const int                      num_targ_vecs
00215   ,VectorBase<Scalar>*const       targ_vecs[]
00216   ,RTOpPack::ReductTarget         *reduct_obj
00217   ,const Index                    first_ele_offset
00218   ,const Index                    sub_dim
00219   ,const Index                    global_offset
00220   ) const
00221 {
00222   applyOp(
00223     NULL
00224     ,op,num_vecs,vecs,num_targ_vecs,targ_vecs,reduct_obj
00225     ,first_ele_offset,sub_dim,global_offset
00226     );
00227 }
00228 
00229 template<class Scalar>
00230 void SpmdVectorBase<Scalar>::acquireDetachedView( const Range1D& rng_in, RTOpPack::ConstSubVectorView<Scalar>* sub_vec ) const
00231 {
00232   if( rng_in == Range1D::Invalid ) {
00233     // Just return an null view
00234     sub_vec->initialize(
00235       rng_in.lbound()    // globalOffset
00236       ,0                 // subDim
00237       ,0                 // values
00238       ,1                 // stride
00239       );
00240     return;
00241   }
00242   const Range1D rng = validateRange(rng_in);
00243   if( rng.lbound() < localOffset_ || localOffset_+localSubDim_-1 < rng.ubound() ) {
00244     // rng consists of off-processor elements so use the default implementation!
00245     VectorDefaultBase<Scalar>::acquireDetachedView(rng_in,sub_vec);
00246     return;
00247   }
00248   // rng consists of all local data so get it!
00249   const Scalar *localValues = NULL;
00250   Index stride = 0;
00251   this->getLocalData(&localValues,&stride);
00252   sub_vec->initialize(
00253     rng.lbound()                               // globalOffset
00254     ,rng.size()                                // subDim
00255     ,localValues+(rng.lbound()-localOffset_)   // values
00256     ,stride                                    // stride
00257     );
00258 }
00259 
00260 template<class Scalar>
00261 void SpmdVectorBase<Scalar>::releaseDetachedView( RTOpPack::ConstSubVectorView<Scalar>* sub_vec ) const
00262 {
00263 #ifdef TEUCHOS_DEBUG
00264   TEST_FOR_EXCEPTION(
00265     sub_vec==NULL || sub_vec->globalOffset() < 0 || sub_vec->globalOffset() + sub_vec->subDim() > globalDim_
00266     ,std::logic_error
00267     ,"SpmdVectorBase<Scalar>::releaseDetachedView(...) : Error, this sub vector was not gotten from acquireDetachedView(...)!"
00268     );
00269 #endif
00270   if( sub_vec->globalOffset() < localOffset_ || localOffset_+localSubDim_ < sub_vec->globalOffset()+sub_vec->subDim() ) {
00271     // Let the default implementation handle it!
00272     VectorDefaultBase<Scalar>::releaseDetachedView(sub_vec);
00273     return;
00274   }
00275   // Nothing to deallocate!
00276   sub_vec->set_uninitialized();
00277 }
00278 
00279 template<class Scalar>
00280 void SpmdVectorBase<Scalar>::acquireDetachedView( const Range1D& rng_in, RTOpPack::SubVectorView<Scalar>* sub_vec )
00281 {
00282   if( rng_in == Range1D::Invalid ) {
00283     // Just return an null view
00284     sub_vec->initialize(
00285       rng_in.lbound()    // globalOffset
00286       ,0                 // subDim
00287       ,0                 // values
00288       ,1                 // stride
00289       );
00290     return;
00291   }
00292   const Range1D rng = validateRange(rng_in);
00293   if( rng.lbound() < localOffset_ || localOffset_+localSubDim_-1 < rng.ubound() ) {
00294     // rng consists of off-processor elements so use the default implementation!
00295     VectorDefaultBase<Scalar>::acquireDetachedView(rng_in,sub_vec);
00296     return;
00297   }
00298   // rng consists of all local data so get it!
00299   Scalar *localValues = NULL;
00300   Index stride = 0;
00301   this->getLocalData(&localValues,&stride);
00302   sub_vec->initialize(
00303     rng.lbound()                               // globalOffset
00304     ,rng.size()                                // subDim
00305     ,localValues+(rng.lbound()-localOffset_)   // values
00306     ,stride                                    // stride
00307     );
00308 }
00309 
00310 template<class Scalar>
00311 void SpmdVectorBase<Scalar>::commitDetachedView( RTOpPack::SubVectorView<Scalar>* sub_vec )
00312 {
00313 #ifdef TEUCHOS_DEBUG
00314   TEST_FOR_EXCEPTION(
00315     sub_vec==NULL || sub_vec->globalOffset() < 0 || sub_vec->globalOffset() + sub_vec->subDim() > globalDim_
00316     ,std::logic_error
00317     ,"SpmdVectorBase<Scalar>::commitDetachedView(...) : Error, this sub vector was not gotten from acquireDetachedView(...)!"
00318     );
00319 #endif
00320   if( sub_vec->globalOffset() < localOffset_ || localOffset_+localSubDim_ < sub_vec->globalOffset()+sub_vec->subDim() ) {
00321     // Let the default implementation handle it!
00322     VectorDefaultBase<Scalar>::commitDetachedView(sub_vec);
00323     return;
00324   }
00325   sub_vec->set_uninitialized();  // Nothing to deallocate!
00326 }
00327 
00328 // protected
00329 
00330 template<class Scalar>
00331 void SpmdVectorBase<Scalar>::updateSpmdSpace()
00332 {
00333   if(globalDim_ == 0) {
00334     const SpmdVectorSpaceBase<Scalar> *spmdSpace = this->spmdSpace().get();
00335     if(spmdSpace) {
00336       globalDim_    = spmdSpace->dim();
00337       localOffset_  = spmdSpace->localOffset();
00338       localSubDim_  = spmdSpace->localSubDim();
00339     }
00340     else {
00341       globalDim_    = 0;
00342       localOffset_  = -1;
00343       localSubDim_  = 0;
00344     }
00345   }
00346 }
00347 
00348 // private
00349 
00350 template<class Scalar>
00351 Range1D SpmdVectorBase<Scalar>::validateRange( const Range1D &rng_in ) const
00352 {
00353   const Range1D rng = Teuchos::full_range(rng_in,0,globalDim_-1);
00354 #ifdef TEUCHOS_DEBUG
00355   TEST_FOR_EXCEPTION(
00356     !(0 <= rng.lbound() && rng.ubound() < globalDim_), std::invalid_argument
00357     ,"SpmdVectorBase<Scalar>::validateRange(...): Error, the range ["
00358     <<rng.lbound()<<","<<rng.ubound()<<"] is not "
00359     "in the range [0,"<<(globalDim_-1)<<"]!"
00360     );
00361 #endif
00362   return rng;
00363 }
00364 
00365 #ifdef THYRA_SPMD_VECTOR_BASE_DUMP
00366 template<class Scalar>
00367 bool SpmdVectorBase<Scalar>::show_dump = false;
00368 #endif // THYRA_SPMD_VECTOR_BASE_DUMP
00369 
00370 } // end namespace Thyra
00371 
00372 #endif // THYRA_SPMD_VECTOR_BASE_HPP

Generated on Thu Sep 18 12:33:03 2008 for Thyra Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1