Kokkos Node API and Local Linear Algebra Kernels Version of the Day
Tsqr_DistTsqr.hpp
Go to the documentation of this file.
00001 //@HEADER
00002 // ************************************************************************
00003 // 
00004 //          Kokkos: Node API and Parallel Node Kernels
00005 //              Copyright (2008) Sandia Corporation
00006 // 
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the U.S. Government retains certain rights in this software.
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 Michael A. Heroux (maherou@sandia.gov) 
00038 // 
00039 // ************************************************************************
00040 //@HEADER
00041 
00045 #ifndef __TSQR_Tsqr_DistTsqr_hpp
00046 #define __TSQR_Tsqr_DistTsqr_hpp
00047 
00048 #include <Tsqr_DistTsqrHelper.hpp>
00049 #include <Tsqr_DistTsqrRB.hpp>
00050 #include <Teuchos_ParameterList.hpp>
00051 #include <Teuchos_ParameterListAcceptorDefaultBase.hpp>
00052 #include <Teuchos_ScalarTraits.hpp>
00053 #include <utility> // std::pair
00054 
00055 
00056 namespace TSQR {
00072   template<class LocalOrdinal, class Scalar>
00073   class DistTsqr : public Teuchos::ParameterListAcceptorDefaultBase {
00074   public:
00075     typedef Scalar scalar_type;
00076     typedef LocalOrdinal ordinal_type;
00077     typedef MatView<ordinal_type, scalar_type > matview_type;
00078     typedef std::vector<std::vector<scalar_type> > VecVec;
00079     typedef std::pair<VecVec, VecVec> FactorOutput;
00080     typedef int rank_type;
00081 
00082   private:
00083     typedef Teuchos::ScalarTraits<Scalar> STS;
00084 
00085   public:
00086 
00093     DistTsqr (const Teuchos::RCP<Teuchos::ParameterList>& plist)
00094     {
00095       setParameterList (plist);
00096     }
00097 
00099     DistTsqr ()
00100     {
00101       setParameterList (Teuchos::null);
00102     }
00103 
00104     void 
00105     setParameterList (const Teuchos::RCP<Teuchos::ParameterList>& plist)
00106     {
00107       using Teuchos::ParameterList;
00108       using Teuchos::parameterList;
00109       using Teuchos::RCP;
00110 
00111       RCP<ParameterList> params = plist.is_null() ? 
00112   parameterList (*getValidParameters()) : plist;
00113 
00114       // Do nothing for now, other than store the list.
00115       this->setMyParamList (params);
00116     }
00117 
00118     Teuchos::RCP<const Teuchos::ParameterList>
00119     getValidParameters() const
00120     {
00121       return Teuchos::parameterList ("DistTsqr"); // Empty list for now.
00122     }
00123     
00128     void init (const Teuchos::RCP<MessengerBase<scalar_type> >& messenger)
00129     {
00130       messenger_ = messenger;
00131       reduceBroadcastImpl_ = 
00132   Teuchos::rcp (new DistTsqrRB<ordinal_type, scalar_type> (messenger_));
00133     }
00134 
00140     rank_type rank() const { 
00141       TEUCHOS_TEST_FOR_EXCEPTION(! ready(), std::logic_error,
00142          "Before using DistTsqr computational methods, "
00143          "you must first call init() with a valid "
00144          "MessengerBase instance.");
00145       return messenger_->rank(); 
00146     }
00147 
00153     rank_type size() const { 
00154       TEUCHOS_TEST_FOR_EXCEPTION(! ready(), std::logic_error,
00155          "Before using DistTsqr computational methods, "
00156          "you must first call init() with a valid "
00157          "MessengerBase instance.");
00158       return messenger_->size(); 
00159     }
00160 
00165     virtual ~DistTsqr () {}
00166 
00175     bool QR_produces_R_factor_with_nonnegative_diagonal () const {
00176       TEUCHOS_TEST_FOR_EXCEPTION(! ready(), std::logic_error,
00177          "Before using DistTsqr computational methods, "
00178          "you must first call init() with a valid "
00179          "MessengerBase instance.");
00180       typedef Combine<ordinal_type, scalar_type> combine_type;
00181       return combine_type::QR_produces_R_factor_with_nonnegative_diagonal() &&
00182   reduceBroadcastImpl_->QR_produces_R_factor_with_nonnegative_diagonal();
00183     }
00184 
00209     void
00210     factorExplicit (matview_type R_mine, 
00211         matview_type Q_mine,
00212         const bool forceNonnegativeDiagonal=false)
00213     {
00214       TEUCHOS_TEST_FOR_EXCEPTION(! ready(), std::logic_error,
00215          "Before using DistTsqr computational methods, "
00216          "you must first call init() with a valid "
00217          "MessengerBase instance.");
00218       reduceBroadcastImpl_->factorExplicit (R_mine, Q_mine, 
00219               forceNonnegativeDiagonal);
00220     }
00221 
00227     void 
00228     getFactorExplicitTimings (std::vector<TimeStats>& stats) const
00229     {
00230       TEUCHOS_TEST_FOR_EXCEPTION(! ready(), std::logic_error,
00231          "Before using DistTsqr computational methods, "
00232          "you must first call init() with a valid "
00233          "MessengerBase instance.");
00234       reduceBroadcastImpl_->getStats (stats);
00235     }
00236 
00242     void
00243     getFactorExplicitTimingLabels (std::vector<std::string>& labels) const
00244     {
00245       TEUCHOS_TEST_FOR_EXCEPTION(! ready(), std::logic_error,
00246          "Before using DistTsqr computational methods, "
00247          "you must first call init() with a valid "
00248          "MessengerBase instance.");
00249       reduceBroadcastImpl_->getStatsLabels (labels);
00250     }
00251 
00275     FactorOutput
00276     factor (matview_type R_mine)
00277     {
00278       TEUCHOS_TEST_FOR_EXCEPTION(! ready(), std::logic_error,
00279          "Before using DistTsqr computational methods, "
00280          "you must first call init() with a valid "
00281          "MessengerBase instance.");
00282       VecVec Q_factors, tau_arrays;
00283       DistTsqrHelper<ordinal_type, scalar_type> helper;
00284       const ordinal_type ncols = R_mine.ncols();
00285 
00286       std::vector< scalar_type > R_local (ncols*ncols);
00287       copy_matrix (ncols, ncols, &R_local[0], ncols, R_mine.get(), R_mine.lda());
00288 
00289       const int P = messenger_->size();
00290       const int my_rank = messenger_->rank();
00291       const int first_tag = 0;
00292       std::vector<scalar_type> work (ncols);
00293       helper.factor_helper (ncols, R_local, my_rank, 0, P-1, first_tag, 
00294           messenger_.get(), Q_factors, tau_arrays, work);
00295       copy_matrix (ncols, ncols, R_mine.get(), R_mine.lda(), &R_local[0], ncols);
00296       return std::make_pair (Q_factors, tau_arrays);
00297     }
00298 
00300     void
00301     apply (const ApplyType& apply_type,
00302      const ordinal_type ncols_C,
00303      const ordinal_type ncols_Q,
00304      scalar_type C_mine[],
00305      const ordinal_type ldc_mine,
00306      const FactorOutput& factor_output)
00307     {
00308       TEUCHOS_TEST_FOR_EXCEPTION(! ready(), std::logic_error,
00309          "Before using DistTsqr computational methods, "
00310          "you must first call init() with a valid "
00311          "MessengerBase instance.");
00312       const bool transposed = apply_type.transposed();
00313       TEUCHOS_TEST_FOR_EXCEPTION(transposed, std::logic_error,
00314          "DistTsqr: Applying Q^T or Q^H has not yet "
00315          "been implemented.");
00316       const int P = messenger_->size();
00317       const int my_rank = messenger_->rank();
00318       const int first_tag = 0;
00319       std::vector<scalar_type> C_other (ncols_C * ncols_C);
00320       std::vector<scalar_type> work (ncols_C);
00321   
00322       const VecVec& Q_factors = factor_output.first;
00323       const VecVec& tau_arrays = factor_output.second;
00324 
00325       // assert (Q_factors.size() == tau_arrays.size());
00326       const int cur_pos = Q_factors.size() - 1;
00327       DistTsqrHelper<ordinal_type, scalar_type> helper;
00328       helper.apply_helper (apply_type, ncols_C, ncols_Q, C_mine, ldc_mine, 
00329          &C_other[0], my_rank, 0, P-1, first_tag, 
00330          messenger_.get(), Q_factors, tau_arrays, cur_pos, 
00331          work);
00332     }
00333 
00335     void
00336     explicit_Q (const ordinal_type ncols_Q,
00337     scalar_type Q_mine[],
00338     const ordinal_type ldq_mine,
00339     const FactorOutput& factor_output)
00340     {
00341       TEUCHOS_TEST_FOR_EXCEPTION(! ready(), std::logic_error,
00342          "Before using DistTsqr computational methods, "
00343          "you must first call init() with a valid "
00344          "MessengerBase instance.");
00345       const int myRank = messenger_->rank ();
00346       fill_matrix (ncols_Q, ncols_Q, Q_mine, ldq_mine, STS::zero());
00347       if (myRank == 0) {
00348   for (ordinal_type j = 0; j < ncols_Q; ++j)
00349     Q_mine[j + j*ldq_mine] = STS::one();
00350       }
00351       apply (ApplyType::NoTranspose, ncols_Q, ncols_Q, 
00352        Q_mine, ldq_mine, factor_output);
00353     }
00354 
00355   private:
00356     Teuchos::RCP<MessengerBase<scalar_type> > messenger_;
00357     Teuchos::RCP<DistTsqrRB<ordinal_type, scalar_type> > reduceBroadcastImpl_;
00358 
00362     bool ready() const {
00363       return ! messenger_.is_null() && ! reduceBroadcastImpl_.is_null();
00364     }
00365   };
00366 
00367 } // namespace TSQR
00368 
00369 #endif // __TSQR_Tsqr_DistTsqr_hpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends