1 // @HEADER 2 // ************************************************************************ 3 // 4 // Rapid Optimization Library (ROL) Package 5 // Copyright (2014) Sandia Corporation 6 // 7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive 8 // license for use of this work by or on behalf of the U.S. Government. 9 // 10 // Redistribution and use in source and binary forms, with or without 11 // modification, are permitted provided that the following conditions are 12 // met: 13 // 14 // 1. Redistributions of source code must retain the above copyright 15 // notice, this list of conditions and the following disclaimer. 16 // 17 // 2. Redistributions in binary form must reproduce the above copyright 18 // notice, this list of conditions and the following disclaimer in the 19 // documentation and/or other materials provided with the distribution. 20 // 21 // 3. Neither the name of the Corporation nor the names of the 22 // contributors may be used to endorse or promote products derived from 23 // this software without specific prior written permission. 24 // 25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY 26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR 28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE 29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, 30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, 31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF 33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING 34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 36 // 37 // Questions? Contact lead developers: 38 // Drew Kouri (dpkouri@sandia.gov) and 39 // Denis Ridzal (dridzal@sandia.gov) 40 // 41 // ************************************************************************ 42 // @HEADER 43 44 #ifndef ROL_PINTVECTORCOMMUNICATION_TPETRA_H 45 #define ROL_PINTVECTORCOMMUNICATION_TPETRA_H 46 47 #include "mpi.h" 48 49 #include "ROL_PinTVectorCommunication.hpp" 50 #include "ROL_StdVector.hpp" 51 #include "ROL_TpetraMultiVector.hpp" 52 53 /** @ingroup la_group 54 \class ROL::PinTVector 55 \brief Defines a parallel in time vector. 56 */ 57 58 namespace ROL { 59 60 template <class Real> 61 class PinTVectorCommunication_Tpetra : public PinTVectorCommunication<Real> { 62 public: 63 /** 64 * \brief Send a vector to a neighboring processor. 65 */ send(MPI_Comm comm,int rank,Vector<Real> & source,int tag=0) const66 void send(MPI_Comm comm,int rank,Vector<Real> & source,int tag=0) const override 67 { 68 auto & tp_source = *dynamic_cast<TpetraMultiVector<Real>&>(source).getVector(); 69 tp_source.sync_host(); 70 71 auto view = tp_source.getLocalViewHost(); 72 73 // int myRank = -1; 74 // MPI_Comm_rank(comm, &myRank); 75 76 MPI_Send(const_cast<Real*>(&view(0,0)),int(view.extent(0)*view.extent(1)),MPI_DOUBLE,rank,tag,comm); 77 } 78 79 /** 80 * \brief Recieve a vector from a neighboring processor. Uses blocking communication. 81 */ recv(MPI_Comm comm,int rank,Vector<Real> & dest,int tag=0) const82 void recv(MPI_Comm comm,int rank,Vector<Real> & dest,int tag=0) const override 83 { 84 auto & tp_dest = *dynamic_cast<TpetraMultiVector<Real>&>(dest).getVector(); 85 auto view = tp_dest.getLocalViewHost(); 86 87 // int myRank = -1; 88 // MPI_Comm_rank(comm, &myRank); 89 90 MPI_Recv(&view(0,0),int(view.extent(0)*view.extent(1)),MPI_DOUBLE,rank,tag,comm,MPI_STATUS_IGNORE); 91 92 // tp_source.template sync<Kokkos::DeviceSpace>(); 93 } 94 95 /** 96 * \brief Recieve a vector from a neighboring processor. Uses blocking communication. 97 */ recvSumInto(MPI_Comm comm,int rank,Vector<Real> & dest,int tag=0) const98 void recvSumInto(MPI_Comm comm,int rank,Vector<Real> & dest,int tag=0) const override 99 { 100 auto & tp_dest = *dynamic_cast<TpetraMultiVector<Real>&>(dest).getVector(); 101 auto view = tp_dest.getLocalViewHost(); 102 103 int myRank = -1; 104 MPI_Comm_rank(comm, &myRank); 105 106 assert(view.extent(1)==1); 107 108 std::vector<Real> buffer(view.extent(0)*view.extent(1),0.0); 109 MPI_Recv(&buffer[0],int(buffer.size()),MPI_DOUBLE,rank,tag,comm,MPI_STATUS_IGNORE); 110 111 for(size_t i=0;i<buffer.size();i++) 112 view(i,0) += buffer[i]; 113 } 114 }; 115 116 } // namespace ROL 117 118 #endif 119