Kokkos Node API and Local Linear Algebra Kernels Version of the Day
Tsqr_printGlobalMatrix.hpp
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 
00042 #ifndef __Tsqr_printGlobalMatrix_hpp
00043 #define __Tsqr_printGlobalMatrix_hpp
00044 
00045 #include <Tsqr_MessengerBase.hpp>
00046 #include <Tsqr_Util.hpp>
00047 #include <Tsqr_Matrix.hpp>
00048 #include <Teuchos_ScalarTraits.hpp>
00049 
00050 #include <limits>
00051 #include <ostream>
00052 #include <stdexcept>
00053 
00054 
00055 namespace TSQR {
00056 
00071   template<class ConstMatrixViewType>
00072   void
00073   printGlobalMatrix (std::ostream& out,
00074          const ConstMatrixViewType& A_local,
00075          MessengerBase<typename ConstMatrixViewType::scalar_type>* const scalarComm,
00076          MessengerBase<typename ConstMatrixViewType::ordinal_type>* const ordinalComm)
00077     {
00078       typedef typename ConstMatrixViewType::ordinal_type LocalOrdinal;
00079       typedef typename ConstMatrixViewType::scalar_type Scalar;
00080       typedef Teuchos::ScalarTraits<Scalar> STS;
00081       using std::endl;
00082 
00083       const int myRank = scalarComm->rank ();
00084       const int nprocs = scalarComm->size ();
00085       const LocalOrdinal nrowsLocal = A_local.nrows();
00086       const LocalOrdinal ncols = A_local.ncols();
00087       const Scalar quiet_NaN = STS::nan();
00088 
00089       if (myRank == 0)
00090   {
00091     // Print the remote matrix data
00092     // out << "Processor " << my_rank << ":" << endl;
00093     print_local_matrix (out, A_local.nrows(), A_local.ncols(), 
00094             A_local.get(), A_local.lda());
00095 
00096     // Space for remote matrix data.  Other processors are allowed
00097     // to have different nrows_local values; we make space as
00098     // necessary.
00099     Matrix<LocalOrdinal, Scalar> A_remote (nrowsLocal, ncols, quiet_NaN);
00100 
00101     // Loop through all the other processors in order.
00102     // Fetch their matrix data and print it.
00103     for (int srcProc = 1; srcProc < nprocs; ++srcProc)
00104       {
00105         // Get processor proc's local matrix dimensions
00106         LocalOrdinal dims[2];
00107         ordinalComm->recv (&dims[0], 2, srcProc, 0);
00108 
00109         // Make space for the remote matrix data.
00110         //
00111         // mfh 13 Oct 2010: Teuchos::OrdinalTraits does not
00112         // currently have this feature.  It's OK to use
00113         // std::numeric_limits, since ordinal types in Trilinos
00114         // are intended to be built-in types (like int or long
00115         // long int).  std::numeric_limits only promises to work
00116         // for built-in types, unless someone has defined an
00117         // appropriate specialization.  Teuchos::ScalarTraits,
00118         // in contrast, has to work for non-built-in Scalar
00119         // types, like ARPREC or QD floating-point numbers.
00120         if (std::numeric_limits<LocalOrdinal>::is_signed)
00121     {
00122       if (dims[0] <= 0 || dims[1] <= 0)
00123         throw std::runtime_error ("Invalid dimensions of remote matrix");
00124     }
00125         else
00126     {
00127       if (dims[0] == 0 || dims[1] == 0)
00128         throw std::runtime_error ("Invalid dimensions of remote matrix");
00129     }
00130         A_remote.reshape (dims[0], dims[1]);
00131 
00132         // Receive the remote matrix data, which we assume is
00133         // stored contiguously.
00134         scalarComm->recv (A_remote.get(), dims[0]*dims[1], srcProc, 0);
00135 
00136         // Print the remote matrix data
00137         // out << "Processor " << proc << ":" << endl;
00138         print_local_matrix (out, dims[0], dims[0], A_remote.get(), A_remote.lda());
00139       }
00140   }
00141       else
00142   {
00143     // Send my local matrix dimensions to proc 0.
00144     int rootProc = 0;
00145     LocalOrdinal dims[2];
00146 
00147     dims[0] = nrowsLocal;
00148     dims[1] = ncols;
00149     ordinalComm->send (dims, 2, rootProc, 0);
00150 
00151     // Create a (contiguous) buffer and copy the data into it.
00152     Matrix< LocalOrdinal, Scalar > A_buf (nrowsLocal, ncols, quiet_NaN);
00153     A_buf.copy (A_local);
00154 
00155     // Send the actual data to proc 0.
00156     scalarComm->send (A_buf.get(), nrowsLocal*ncols, rootProc, 0);
00157   }
00158       scalarComm->barrier ();
00159     }
00160 
00161 } // namespace TSQR
00162 
00163 #endif // __Tsqr_printGlobalMatrix_hpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends