Kokkos Node API and Local Linear Algebra Kernels Version of the Day
CrsMatrix_DefaultMultiply.cpp
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 #include <Teuchos_UnitTestHarness.hpp>
00043 #include <Teuchos_TimeMonitor.hpp>
00044 #include <Teuchos_Time.hpp>
00045 #include <Teuchos_TypeNameTraits.hpp>
00046 #include <Teuchos_ScalarTraits.hpp>
00047 
00048 #include "Kokkos_ConfigDefs.hpp"
00049 #include "Kokkos_DefaultArithmetic.hpp"
00050 #include "Kokkos_DefaultKernels.hpp"
00051 #include "Kokkos_Version.hpp"
00052 
00053 #include "Kokkos_SerialNode.hpp"
00054 #ifdef HAVE_KOKKOS_TBB
00055 #include "Kokkos_TBBNode.hpp"
00056 #endif
00057 #ifdef HAVE_KOKKOS_THREADPOOL
00058 #include "Kokkos_TPINode.hpp"
00059 #endif
00060 #ifdef HAVE_KOKKOS_THRUST
00061 #include "Kokkos_ThrustGPUNode.hpp"
00062 #endif
00063 
00064 namespace {
00065 
00066   using Kokkos::MultiVector;
00067   using Kokkos::CrsMatrix;
00068   using Kokkos::CrsGraph;
00069   using Kokkos::DefaultArithmetic;
00070   using Kokkos::DefaultKernels;
00071   using Kokkos::SerialNode;
00072   using Teuchos::ArrayRCP;
00073   using Teuchos::RCP;
00074   using Teuchos::rcp;
00075   using Teuchos::null;
00076   using std::endl;
00077 
00078   RCP<SerialNode> snode;
00079 #ifdef HAVE_KOKKOS_TBB
00080   using Kokkos::TBBNode;
00081   RCP<TBBNode> tbbnode;
00082 #endif
00083 #ifdef HAVE_KOKKOS_THREADPOOL
00084   using Kokkos::TPINode;
00085   RCP<TPINode> tpinode;
00086 #endif
00087 #ifdef HAVE_KOKKOS_THRUST
00088   using Kokkos::ThrustGPUNode;
00089   RCP<ThrustGPUNode> thrustnode;
00090 #endif
00091 
00092   int N = 1000;
00093 
00094   TEUCHOS_STATIC_SETUP()
00095   {
00096     Teuchos::CommandLineProcessor &clp = Teuchos::UnitTestRepository::getCLP();
00097     clp.addOutputSetupOptions(true);
00098     clp.setOption("test-size",&N,"Vector length for tests.");
00099   }
00100 
00101   template <class Node>
00102   RCP<Node> getNode() {
00103     assert(false);
00104   }
00105 
00106   template <>
00107   RCP<SerialNode> getNode<SerialNode>() {
00108     if (snode == null) {
00109       Teuchos::ParameterList pl;
00110       snode = rcp(new SerialNode(pl));
00111     }
00112     return snode;
00113   }
00114 
00115 #ifdef HAVE_KOKKOS_TBB
00116   template <>
00117   RCP<TBBNode> getNode<TBBNode>() {
00118     if (tbbnode == null) {
00119       Teuchos::ParameterList pl;
00120       pl.set<int>("Num Threads",0);
00121       tbbnode = rcp(new TBBNode(pl));
00122     }
00123     return tbbnode;
00124   }
00125 #endif
00126 
00127 #ifdef HAVE_KOKKOS_THREADPOOL
00128   template <>
00129   RCP<TPINode> getNode<TPINode>() {
00130     if (tpinode == null) {
00131       Teuchos::ParameterList pl;
00132       pl.set<int>("Num Threads",0);
00133       tpinode = rcp(new TPINode(pl));
00134     }
00135     return tpinode;
00136   }
00137 #endif
00138 
00139 #ifdef HAVE_KOKKOS_THRUST
00140   template <>
00141   RCP<ThrustGPUNode> getNode<ThrustGPUNode>() {
00142     if (thrustnode == null) {
00143       Teuchos::ParameterList pl;
00144       pl.set<int>("Num Threads",0);
00145       pl.set<int>("Verbose",1);
00146       thrustnode = rcp(new ThrustGPUNode(pl));
00147     }
00148     return thrustnode;
00149   }
00150 #endif
00151 
00152   //
00153   // UNIT TESTS
00154   // 
00155 
00156   TEUCHOS_UNIT_TEST_TEMPLATE_3_DECL( CrsMatrix, SparseMultiply, Ordinal, Scalar, Node )
00157   {
00158     RCP<Node> node = getNode<Node>();
00159     typedef typename DefaultKernels<Scalar,Ordinal,Node>::SparseOps DSM;
00160     typedef CrsGraph<Ordinal,Node,DSM>                             GRPH;
00161     typedef CrsMatrix<Scalar,Ordinal,Node,DSM>                      MAT;
00162     typedef MultiVector<Scalar,Node>                                 MV;
00163     typedef Teuchos::ScalarTraits<Scalar>                            ST;
00164     const Scalar ONE = ST::one(),
00165                 ZERO = ST::zero();
00166     // generate tridiagonal matrix:
00167     // [ 2 -1                   ]
00168     // [-1  3  -1               ]
00169     // [   -1   3  -1           ]
00170     // [                        ]
00171     // [                -1  3 -1]
00172     // [                   -1  2]
00173     if (N<2) return;
00174     GRPH G(N,node);
00175     MAT  A(G);
00176     // allocate buffers for offsets, indices and values
00177     const size_t totalNNZ = 3*N - 2;
00178     ArrayRCP<size_t> offsets(N+1);
00179     ArrayRCP<Ordinal>   inds(totalNNZ);
00180     ArrayRCP<Scalar>    vals(totalNNZ);
00181     // fill the buffers on the host
00182     {
00183       size_t NNZsofar = 0;
00184       offsets[0] = NNZsofar;
00185       inds[NNZsofar] = 0; inds[NNZsofar+1] =  1;
00186       vals[NNZsofar] = 2; vals[NNZsofar+1] = -1;
00187       NNZsofar += 2;
00188       for (int i=1; i != N-1; ++i) {
00189         offsets[i] = NNZsofar;
00190         inds[NNZsofar] = i-1; inds[NNZsofar+1] = i; inds[NNZsofar+2] = i+1;
00191         vals[NNZsofar] =  -1; vals[NNZsofar+1] = 3; vals[NNZsofar+2] =  -1;
00192         NNZsofar += 3;
00193       }
00194       offsets[N-1] = NNZsofar;
00195       inds[NNZsofar] = N-2; inds[NNZsofar+1] = N-1;
00196       vals[NNZsofar] =  -1; vals[NNZsofar+1] = 2;
00197       NNZsofar += 2;
00198       offsets[N]   = NNZsofar;
00199       TEUCHOS_TEST_FOR_EXCEPT(NNZsofar != totalNNZ);
00200     }
00201     G.set1DStructure(inds, offsets, offsets.persistingView(1,N));
00202     offsets = Teuchos::null;
00203     inds    = Teuchos::null;
00204     A.set1DValues(vals);
00205     vals    = Teuchos::null;
00206     A.finalize(true);
00207     typename DSM::template rebind<Scalar>::other dsm(node);
00208     out << "Testing with sparse ops: " << Teuchos::typeName(dsm) << std::endl;
00209     dsm.initializeStructure(G);
00210     dsm.initializeValues(A);
00211 
00212     ArrayRCP<Scalar> xdat, axdat;
00213     xdat  = node->template allocBuffer<Scalar>(N);
00214     axdat = node->template allocBuffer<Scalar>(N);
00215     MV X(node), AX(node);
00216     X.initializeValues( N,1, xdat,N);
00217     AX.initializeValues(N,1,axdat,N);
00218     DefaultArithmetic<MV>::Init( X,1);
00219     dsm.multiply(Teuchos::NO_TRANS,ONE,X,AX);
00220     // AX should be all ones
00221     {
00222       ArrayRCP<const Scalar> axview = node->template viewBuffer<Scalar>(N,axdat);
00223       Scalar err = ZERO;
00224       for (int i=0; i<N; ++i) {
00225         err += ST::magnitude(ONE - axview[i]);
00226       }
00227       TEST_EQUALITY_CONST(err, ZERO);
00228     }
00229     // do that same multiplication, testing alpha=-1 and beta=1 accumulation
00230     dsm.multiply(Teuchos::NO_TRANS,-ONE,X,ONE,AX);
00231     // AX should be zero
00232     {
00233       ArrayRCP<const Scalar> axview = node->template viewBuffer<Scalar>(N,axdat);
00234       Scalar err = ZERO;
00235       for (int i=0; i<N; ++i) {
00236         err += ST::magnitude(axview[i]);
00237       }
00238       TEST_EQUALITY_CONST(err, ZERO);
00239     }
00240     xdat = null;
00241     axdat = null;
00242   }
00243 
00244 #include "CrsMatrix_DefaultMultiplyTests.hpp"
00245 
00246 #define ALL_UNIT_TESTS_ORDINAL_SCALAR_NODE( ORDINAL, SCALAR, NODE ) \
00247       TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( CrsMatrix,        SparseMultiply, ORDINAL, SCALAR, NODE ) \
00248       TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( DefaultSparseOps, NodeTest,       ORDINAL, SCALAR, NODE ) \
00249       TEUCHOS_UNIT_TEST_TEMPLATE_3_INSTANT( DefaultSparseOps, ResubmitMatrix, ORDINAL, SCALAR, NODE )
00250 
00251 #define UNIT_TEST_SERIALNODE(ORDINAL, SCALAR) \
00252       ALL_UNIT_TESTS_ORDINAL_SCALAR_NODE( ORDINAL, SCALAR, SerialNode )
00253 
00254 #ifdef HAVE_KOKKOS_TBB
00255 #define UNIT_TEST_TBBNODE(ORDINAL, SCALAR) \
00256       ALL_UNIT_TESTS_ORDINAL_SCALAR_NODE( ORDINAL, SCALAR, TBBNode )
00257 #else
00258 #define UNIT_TEST_TBBNODE(ORDINAL, SCALAR)
00259 #endif
00260 
00261 #ifdef HAVE_KOKKOS_THREADPOOL
00262 #define UNIT_TEST_TPINODE(ORDINAL, SCALAR) \
00263       ALL_UNIT_TESTS_ORDINAL_SCALAR_NODE( ORDINAL, SCALAR, TPINode )
00264 #else
00265 #define UNIT_TEST_TPINODE(ORDINAL, SCALAR)
00266 #endif
00267 
00268 #ifdef HAVE_KOKKOS_THRUST
00269 #define UNIT_TEST_THRUSTGPUNODE(ORDINAL, SCALAR) \
00270       ALL_UNIT_TESTS_ORDINAL_SCALAR_NODE( ORDINAL, SCALAR, ThrustGPUNode )
00271 #else
00272 #define UNIT_TEST_THRUSTGPUNODE(ORDINAL, SCALAR)
00273 #endif
00274 
00275 #define UNIT_TEST_GROUP_ORDINAL_SCALAR( ORDINAL, SCALAR ) \
00276         UNIT_TEST_SERIALNODE( ORDINAL, SCALAR ) \
00277         UNIT_TEST_TBBNODE( ORDINAL, SCALAR ) \
00278         UNIT_TEST_TPINODE( ORDINAL, SCALAR ) \
00279         UNIT_TEST_THRUSTGPUNODE( ORDINAL, SCALAR )
00280 
00281 #define UNIT_TEST_GROUP_ORDINAL( ORDINAL ) \
00282         UNIT_TEST_GROUP_ORDINAL_SCALAR(ORDINAL, int) \
00283         UNIT_TEST_GROUP_ORDINAL_SCALAR(ORDINAL, float)
00284 
00285      UNIT_TEST_GROUP_ORDINAL(int)
00286      typedef short int ShortInt; UNIT_TEST_GROUP_ORDINAL(ShortInt)
00287 
00288 }
00289 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends