Thyra Version of the Day
Thyra_DefaultSpmdMultiVector_def.hpp
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 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Roscoe A. Bartlett (bartlettra@ornl.gov) 
00038 // 
00039 // ***********************************************************************
00040 // @HEADER
00041 
00042 #ifndef THYRA_DEFAULT_SPMD_MULTI_VECTOR_DEF_HPP
00043 #define THYRA_DEFAULT_SPMD_MULTI_VECTOR_DEF_HPP
00044 
00045 // Define to make some verbose output
00046 //#define THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
00047 
00048 #include "Thyra_DefaultSpmdMultiVector_decl.hpp"
00049 #include "Thyra_SpmdMultiVectorBase.hpp"
00050 #include "Thyra_VectorSpaceFactoryBase.hpp"
00051 #include "Thyra_DefaultSpmdVector.hpp"
00052 #include "Teuchos_Assert.hpp"
00053 
00054 
00055 namespace Thyra {
00056 
00057 
00058 //
00059 // Simple utility class to copy back multi-vector entries
00060 //
00061 // ToDo: Refactor above to not use raw pointers!
00062 //
00063 template<class Scalar>
00064 class CopyBackSpmdMultiVectorEntries {
00065 public:
00066   CopyBackSpmdMultiVectorEntries(
00067     const ArrayView<const int> &cols,
00068     const ArrayRCP<const Scalar> &localValuesView, const Ordinal localSubDim,
00069     const ArrayRCP<Scalar> &localValues, const Ordinal leadingDim
00070     )
00071     : cols_(cols), localValuesView_(localValuesView), localSubDim_(localSubDim),
00072       localValues_(localValues), leadingDim_(leadingDim)
00073     {}
00074   ~CopyBackSpmdMultiVectorEntries()
00075     {
00076       typedef typename ArrayRCP<const Scalar>::const_iterator const_itr_t;
00077       typedef typename ArrayRCP<Scalar>::iterator itr_t;
00078       // Copy from contiguous storage column by column
00079       if (localValues_.strong_count()) {
00080         const int numCols = cols_.size();
00081         const const_itr_t lvv = localValuesView_.begin();
00082         const itr_t lv = localValues_.begin();
00083         for (int k = 0; k < numCols; ++k) {
00084           const int col_k = cols_[k];
00085           const const_itr_t lvv_k = lvv + localSubDim_*k;
00086           const itr_t lv_k = lv + leadingDim_*col_k;
00087           std::copy( lvv_k, lvv_k + localSubDim_, lv_k );
00088         }
00089       }
00090 #ifdef THYRA_DEBUG
00091       else {
00092         ++DefaultSpmdMultiVector<Scalar>::numSkipCopyBack;
00093       }
00094 #endif // THYRA_DEBUG
00095     }
00096 private:
00097   Array<int> cols_;
00098   ArrayRCP<const Scalar> localValuesView_;
00099   Ordinal localSubDim_;
00100   ArrayRCP<Scalar> localValues_;
00101   Ordinal leadingDim_;
00102   // Not defined and not to be called
00103   CopyBackSpmdMultiVectorEntries();
00104   CopyBackSpmdMultiVectorEntries(const CopyBackSpmdMultiVectorEntries&);
00105   CopyBackSpmdMultiVectorEntries& operator=(const CopyBackSpmdMultiVectorEntries&);
00106 };
00107 
00108 
00109 template<class Scalar>
00110 RCP<CopyBackSpmdMultiVectorEntries<Scalar> >
00111 copyBackSpmdMultiVectorEntries(
00112   const ArrayView<const int> &cols,
00113   const ArrayRCP<const Scalar> &localValuesView, const Ordinal localSubDim,
00114   const ArrayRCP<Scalar> &localValues, const Ordinal leadingDim
00115   )
00116 {
00117   return Teuchos::rcp(
00118     new CopyBackSpmdMultiVectorEntries<Scalar>(
00119       cols, localValuesView, localSubDim, localValues, leadingDim
00120       )
00121     );
00122 }
00123 
00124 
00125 //
00126 // DefaultSpmdMultiVector
00127 //
00128 
00129 
00130 #ifdef THYRA_DEBUG
00131 template<class Scalar>
00132 int DefaultSpmdMultiVector<Scalar>::numSkipCopyBack(0);
00133 #endif
00134 
00135 
00136 // Constructors/initializers/accessors
00137 
00138 
00139 template<class Scalar>
00140 DefaultSpmdMultiVector<Scalar>::DefaultSpmdMultiVector()
00141   :leadingDim_(0)
00142 {}
00143 
00144 
00145 template<class Scalar>
00146 DefaultSpmdMultiVector<Scalar>::DefaultSpmdMultiVector(
00147   const RCP<const SpmdVectorSpaceBase<Scalar> > &spmdRangeSpace
00148   ,const RCP<const ScalarProdVectorSpaceBase<Scalar> > &domainSpace
00149   )
00150 {
00151   initialize(spmdRangeSpace,domainSpace);
00152 }
00153 
00154 
00155 template<class Scalar>
00156 DefaultSpmdMultiVector<Scalar>::DefaultSpmdMultiVector(
00157   const RCP<const SpmdVectorSpaceBase<Scalar> > &spmdRangeSpace
00158   ,const RCP<const ScalarProdVectorSpaceBase<Scalar> > &domainSpace
00159   ,const ArrayRCP<Scalar> &localValues
00160   ,const Ordinal leadingDim
00161   )
00162 {
00163   initialize(spmdRangeSpace,domainSpace,localValues,leadingDim);
00164 }
00165 
00166 
00167 template<class Scalar>
00168 void DefaultSpmdMultiVector<Scalar>::initialize(
00169   const RCP<const SpmdVectorSpaceBase<Scalar> > &spmdRangeSpace
00170   ,const RCP<const ScalarProdVectorSpaceBase<Scalar> > &domainSpace
00171   )
00172 {
00173   const Ordinal localSubDim = spmdRangeSpace->localSubDim();
00174   ArrayRCP<Scalar> values;
00175   if (localSubDim)
00176     values = Teuchos::arcp<Scalar>(localSubDim * domainSpace->dim());
00177   initialize(spmdRangeSpace, domainSpace, values, localSubDim);
00178 }
00179 
00180 
00181 template<class Scalar>
00182 void DefaultSpmdMultiVector<Scalar>::initialize(
00183   const RCP<const SpmdVectorSpaceBase<Scalar> > &spmdRangeSpace,
00184   const RCP<const ScalarProdVectorSpaceBase<Scalar> > &domainSpace,
00185   const ArrayRCP<Scalar> &localValues,
00186   const Ordinal leadingDim_in
00187   )
00188 {
00189   const Ordinal leadingDim =
00190     (leadingDim_in >= 0 ? leadingDim_in : spmdRangeSpace->localSubDim());
00191 #ifdef TEUCHOS_DEBUG
00192   TEUCHOS_ASSERT(!is_null(spmdRangeSpace));
00193   TEUCHOS_ASSERT(!is_null(domainSpace));
00194   if (spmdRangeSpace->dim()) {
00195     TEUCHOS_ASSERT(!is_null(localValues));
00196   }
00197   TEUCHOS_ASSERT_INEQUALITY(leadingDim, >=, spmdRangeSpace->localSubDim());
00198 #endif
00199   spmdRangeSpace_ = spmdRangeSpace;
00200   domainSpace_ = domainSpace;
00201   localValues_ = localValues;
00202   leadingDim_ = leadingDim;
00203   this->updateSpmdSpace();
00204 }
00205 
00206 
00207 template<class Scalar>
00208 void DefaultSpmdMultiVector<Scalar>::uninitialize(
00209   RCP<const SpmdVectorSpaceBase<Scalar> > *spmdRangeSpace
00210   ,RCP<const ScalarProdVectorSpaceBase<Scalar> > *domainSpace
00211   ,ArrayRCP<Scalar> *localValues
00212   ,Ordinal *leadingDim
00213   )
00214 {
00215   if(spmdRangeSpace) *spmdRangeSpace = spmdRangeSpace_;
00216   if(domainSpace) *domainSpace = domainSpace_;
00217   if(localValues) *localValues = localValues_;
00218   if(leadingDim) *leadingDim = leadingDim_;
00219 
00220   spmdRangeSpace_ = Teuchos::null;
00221   domainSpace_ = Teuchos::null;
00222   localValues_ = Teuchos::null;
00223   leadingDim_ = 0;
00224 
00225   this->updateSpmdSpace();
00226 }
00227 
00228 
00229 template<class Scalar>
00230 RCP< const ScalarProdVectorSpaceBase<Scalar> >
00231 DefaultSpmdMultiVector<Scalar>::domainScalarProdVecSpc() const
00232 {
00233 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
00234   std::cerr << "\nSpmdMultiVectorStd<Scalar>::domainScalarProdVecSpc() const called!\n";
00235 #endif
00236   return domainSpace_;
00237 }
00238 
00239 
00240 // Overridden public functions from SpmdMultiVectorBase
00241 
00242 
00243 template<class Scalar>
00244 RCP<const SpmdVectorSpaceBase<Scalar> >
00245 DefaultSpmdMultiVector<Scalar>::spmdSpace() const
00246 {
00247 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
00248   std::cerr << "\nSpmdMultiVectorStd<Scalar>::spmdSpace() const called!\n";
00249 #endif
00250   return spmdRangeSpace_;
00251 }
00252 
00253 
00254 // Overridden protected functions from MultiVectorBase
00255 
00256 
00257 template<class Scalar>
00258 RCP<VectorBase<Scalar> >
00259 DefaultSpmdMultiVector<Scalar>::nonconstColImpl(Ordinal j)
00260 {
00261 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
00262   std::cerr << "\nSpmdMultiVectorStd<Scalar>::col() called!\n";
00263 #endif
00264 #ifdef TEUCHOS_DEBUG
00265   TEUCHOS_TEST_FOR_EXCEPT( !( 0 <= j && j < this->domain()->dim() ) );
00266 #endif
00267   return Teuchos::rcp(
00268     new DefaultSpmdVector<Scalar>(
00269       spmdRangeSpace_,
00270       localValues_.persistingView(j*leadingDim_,spmdRangeSpace_->localSubDim()),
00271       1
00272       )
00273     );
00274   //return Teuchos::rcp(new DefaultVectorMultiVector<Scalar>(subView(Range1D(j,j))));
00275 }
00276 
00277 
00278 template<class Scalar>
00279 RCP<const MultiVectorBase<Scalar> >
00280 DefaultSpmdMultiVector<Scalar>::contigSubViewImpl(
00281   const Range1D& col_rng_in
00282   ) const
00283 {
00284   const Range1D colRng = this->validateColRange(col_rng_in);
00285   return Teuchos::rcp(
00286     new DefaultSpmdMultiVector<Scalar>(
00287       spmdRangeSpace_,
00288       Teuchos::rcp_dynamic_cast<const ScalarProdVectorSpaceBase<Scalar> >(
00289         spmdRangeSpace_->smallVecSpcFcty()->createVecSpc(colRng.size())
00290         ,true
00291         ),
00292       localValues_.persistingView(colRng.lbound()*leadingDim_,colRng.size()*spmdRangeSpace_->localSubDim()),
00293       leadingDim_
00294       )
00295     );
00296 }
00297 
00298 
00299 template<class Scalar>
00300 RCP<MultiVectorBase<Scalar> >
00301 DefaultSpmdMultiVector<Scalar>::nonconstContigSubViewImpl(
00302   const Range1D& col_rng_in
00303   )
00304 {
00305   return Teuchos::rcp_const_cast<MultiVectorBase<Scalar> >(
00306     this->contigSubViewImpl(col_rng_in));
00307   // Have the nonconst version call the const version.  Note that in this case
00308   // we just need to take the const off of the returned MultiVectorBase object
00309   // because the localValues is already handled as nonconst.  This is the
00310   // perfect instance where the advice in Item 3 in "Effective C++ 3rd
00311   // edition" where Scott Meyers recommends having the nonconst version call
00312   // the const version.
00313 }
00314 
00315 
00316 template<class Scalar>
00317 RCP<const MultiVectorBase<Scalar> >
00318 DefaultSpmdMultiVector<Scalar>::nonContigSubViewImpl(
00319   const ArrayView<const int> &cols
00320   ) const
00321 {
00322   THYRA_DEBUG_ASSERT_MV_COLS("nonContigSubViewImpl(cols)", cols);
00323   const int numCols = cols.size();
00324   const ArrayRCP<Scalar> localValuesView = createContiguousCopy(cols);
00325   return defaultSpmdMultiVector<Scalar>(
00326     spmdRangeSpace_,
00327     createSmallScalarProdVectorSpaceBase<Scalar>(spmdRangeSpace_, numCols),
00328     localValuesView
00329     );
00330 }
00331 
00332 
00333 template<class Scalar>
00334 RCP<MultiVectorBase<Scalar> >
00335 DefaultSpmdMultiVector<Scalar>::nonconstNonContigSubViewImpl(
00336   const ArrayView<const int> &cols )
00337 {
00338   THYRA_DEBUG_ASSERT_MV_COLS("nonContigSubViewImpl(cols)", cols);
00339   const int numCols = cols.size();
00340   const ArrayRCP<Scalar> localValuesView = createContiguousCopy(cols);
00341   const Ordinal localSubDim = spmdRangeSpace_->localSubDim();
00342   RCP<CopyBackSpmdMultiVectorEntries<Scalar> > copyBackView =
00343     copyBackSpmdMultiVectorEntries<Scalar>(cols, localValuesView.getConst(),
00344       localSubDim, localValues_.create_weak(), leadingDim_);
00345   return Teuchos::rcpWithEmbeddedObjPreDestroy(
00346     new DefaultSpmdMultiVector<Scalar>(
00347       spmdRangeSpace_,
00348       createSmallScalarProdVectorSpaceBase<Scalar>(spmdRangeSpace_, numCols),
00349       localValuesView),
00350     copyBackView
00351     );
00352 }
00353 
00354 
00355 // Overridden protected members from SpmdMultiVectorBase
00356 
00357 
00358 template<class Scalar>
00359 void DefaultSpmdMultiVector<Scalar>::getNonconstLocalDataImpl(
00360   const Ptr<ArrayRCP<Scalar> > &localValues, const Ptr<Ordinal> &leadingDim
00361   )
00362 {
00363 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
00364   std::cerr << "\nSpmdMultiVectorStd<Scalar>::getLocalDataImpl() called!\n";
00365 #endif
00366   *localValues = localValues_;
00367   *leadingDim = leadingDim_;
00368 }
00369 
00370 
00371 template<class Scalar>
00372 void DefaultSpmdMultiVector<Scalar>::getLocalDataImpl(
00373   const Ptr<ArrayRCP<const Scalar> > &localValues, const Ptr<Ordinal> &leadingDim
00374   ) const
00375 {
00376 #ifdef THYRA_DEFAULT_SPMD_MULTI_VECTOR_VERBOSE_TO_ERROR_OUT
00377   std::cerr << "\nSpmdMultiVectorStd<Scalar>::getLocalData() called!\n";
00378 #endif
00379   *localValues = localValues_;
00380   *leadingDim = leadingDim_;
00381 }
00382 
00383 
00384 // private
00385 
00386 
00387 template<class Scalar>
00388 ArrayRCP<Scalar>
00389 DefaultSpmdMultiVector<Scalar>::createContiguousCopy(
00390   const ArrayView<const int> &cols ) const
00391 {
00392   typedef typename ArrayRCP<Scalar>::const_iterator const_itr_t;
00393   typedef typename ArrayRCP<Scalar>::iterator itr_t;
00394   const int numCols = cols.size();
00395   const Ordinal localSubDim = spmdRangeSpace_->localSubDim();
00396   ArrayRCP<Scalar> localValuesView = Teuchos::arcp<Scalar>(numCols*localSubDim);
00397   // Copy to contiguous storage column by column
00398   const const_itr_t lv = localValues_.begin();
00399   const itr_t lvv = localValuesView.begin();
00400   for (int k = 0; k < numCols; ++k) {
00401     const int col_k = cols[k];
00402     const const_itr_t lv_k = lv + leadingDim_*col_k;
00403     const itr_t lvv_k = lvv + localSubDim*k;
00404     std::copy(lv_k, lv_k+localSubDim, lvv_k);
00405   }
00406   return localValuesView;
00407 }
00408 
00409 
00410 } // end namespace Thyra
00411 
00412 
00413 #endif // THYRA_DEFAULT_SPMD_MULTI_VECTOR_DEF_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines