Hello world

Basic example using Eigen only

Clone pressio-ops and download Eigen-3.4
Copy the snippet below into a main.cpp somewhere:

#include <pressio/ops.hpp>
#include <Eigen/Dense>
#include <iostream>

int main(){
 using namespace pressio;

 Eigen::VectorXd x(3), y(3);
 x << 1.0, 2.0, 3.0;
 y.setConstant(2.0);

 // BLAS-like dot
 const auto d = ops::dot(x, y); // 1*2 + 2*2 + 3*2 = 12

 // y = 1.5*x - y, so y becomes [-0.5, 1.0, 2.5]
 ops::update(y, -1., x, 1.5);

 // z = A * x  (matrix–vector product)
 Eigen::MatrixXd A = Eigen::MatrixXd::Identity(3,3);
 Eigen::VectorXd z(3);
 ops::product(pressio::nontranspose{}, 1., A, x, 0., z);

 std::cout << "dot(x,y) = " << d << '\n';
 std::cout << "y after axpby = " << y.transpose() << '\n';
 std::cout << "z = " << z.transpose() << '\n';

 return 0;
}

Build/run (requires C++17):

g++ -D PRESSIO_ENABLE_TPL_EIGEN -I <path-to-pressio-ops>/include/ -I <path-to-eigen>/eigen-3.4.0 main.cpp
./a.out

For more details about configuring and option, see here.

Basic example using Trilinos

This assumes you have trilinos installed somewhere.

#include <pressio/ops.hpp>
#include <Eigen/Dense>
#include <Tpetra_Core.hpp>
#include <Tpetra_Map.hpp>
#include <Tpetra_Vector.hpp>
#include <Tpetra_MultiVector.hpp>
#include <Teuchos_RCP.hpp>
#include <iostream>

int main(int argc, char** argv)
{
 Tpetra::ScopeGuard scope(&argc, &argv);
 {
   using map_t = Tpetra::Map<>;
   using vec_t = Tpetra::Vector<>;
   using mv_t  = Tpetra::MultiVector<>;

   auto comm = Tpetra::getDefaultComm();
   const int rank = comm->getRank();

   // Global problem size N (distributed across ranks)
   const Tpetra::global_size_t N = 9;

   // Create a contiguous, uniformly distributed map
   using LO = map_t::local_ordinal_type;
   using GO = map_t::global_ordinal_type;
   auto map = Tpetra::createUniformContigMap<LO,GO>(N, comm);

   // vectors
   vec_t x(map), y(map), z(map);
   pressio::ops::fill(x, 1.0);
   pressio::ops::fill(y, 2.0);
   pressio::ops::fill(z, 1.0);

   // Global dot product via pressio-ops
   const auto d = pressio::ops::dot(x, y);

   // MultiVector with 2 columns (row map as x,y,z)
   mv_t MV(map, 2);

   // Obtain column views (each behaves like a Vector in ops)
   auto col0 = pressio::column(MV, 0);
   auto col1 = pressio::column(MV, 1);
   // these modify the underlying MV because are mutable views
   pressio::ops::fill(col0, 3.);
   pressio::ops::fill(col1, 4.);

   // dot between cols
   const auto dcols = pressio::ops::dot(col0, col1);

   // z = z + 0.5 * col1
   pressio::ops::update(z, 1., col1, 0.5);
   auto dotz = pressio::ops::dot(col0, z);

   // perform mat-vec, storing result into Eigen vector
   // v0 = MV^T * z
   Eigen::VectorXd v0(2);
   pressio::ops::fill(v0, 0.);
   pressio::ops::product(pressio::transpose{}, 1., MV, z, 0., v0);

   if (rank == 0){
     std::cout << "dot(x,y) = " << d << "\n";
     std::cout << "dot(col0, col1) = " << dcols << "\n";
     std::cout << "dotz = " << dotz << "\n";
     std::cout << "v0 = " << v0 << "\n";
   }
 }
 return 0;
}