#9: add elementwise_inverse op

include/pressio/ops.hpp

@@ -88,6 +88,7 @@ template<class ...> struct matching_extents;
 #include "ops/eigen/ops_rank1_update.hpp"
 #include "ops/eigen/ops_rank2_update.hpp"
 #include "ops/eigen/ops_elementwise_multiply.hpp"
+#include "ops/eigen/ops_elementwise_inverse.hpp"
 #include "ops/eigen/ops_level2.hpp"
 #include "ops/eigen/ops_level3.hpp"
 #endif

include/pressio/ops/eigen/ops_elementwise_inverse.hpp

@@ -0,0 +1,87 @@
+/*
+//@HEADER
+// ************************************************************************
+//
+// ops_elementwise_inverse.hpp
+//                     		  Pressio
+//                             Copyright 2019
+//    National Technology & Engineering Solutions of Sandia, LLC (NTESS)
+//
+// Under the terms of Contract DE-NA0003525 with NTESS, the
+// U.S. Government retains certain rights in this software.
+//
+// Pressio is licensed under BSD-3-Clause terms of use:
+//
+// Redistribution and use in source and binary forms, with or without
+// modification, are permitted provided that the following conditions
+// are met:
+//
+// 1. Redistributions of source code must retain the above copyright
+// notice, this list of conditions and the following disclaimer.
+//
+// 2. Redistributions in binary form must reproduce the above copyright
+// notice, this list of conditions and the following disclaimer in the
+// documentation and/or other materials provided with the distribution.
+//
+// 3. Neither the name of the copyright holder nor the names of its
+// contributors may be used to endorse or promote products derived
+// from this software without specific prior written permission.
+//
+// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
+// "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
+// LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
+// FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
+// COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
+// INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
+// (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+// SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+// HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
+// STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING
+// IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+// POSSIBILITY OF SUCH DAMAGE.
+//
+// Questions? Contact Francesco Rizzi ([email protected])
+//
+// ************************************************************************
+//@HEADER
+*/
+
+#ifndef PRESSIOOPS_OPS_EIGEN_OPS_ELEMENTWISE_INVERSE_HPP_
+#define PRESSIOOPS_OPS_EIGEN_OPS_ELEMENTWISE_INVERSE_HPP_
+
+namespace pressio{ namespace ops{
+
+//----------------------------------------------------------------------
+// computing elementwise:  y = 1/z
+//----------------------------------------------------------------------
+
+template <class T1, class T2>
+std::enable_if_t<
+  // common elementwise_inverse constraints
+     ::pressio::Traits<T1>::rank == 1
+  && ::pressio::Traits<T2>::rank == 1
+  // TPL/container specific
+  && (::pressio::is_native_container_eigen<T1>::value
+  || ::pressio::is_expression_acting_on_eigen<T1>::value)
+  && (::pressio::is_native_container_eigen<T2>::value
+  || ::pressio::is_expression_acting_on_eigen<T2>::value)
+  // scalar compatibility
+  && ::pressio::all_have_traits_and_same_scalar<T1, T2>::value
+  /* constrained via is_convertible because the impl is using
+     Eigen native operations which are based on expressions and require
+     coefficients to be convertible to scalar types of the vector/matrices */
+   && (std::is_floating_point<typename ::pressio::Traits<T1>::scalar_type>::value
+   || std::is_integral<typename ::pressio::Traits<T1>::scalar_type>::value)
+  >
+elementwise_inverse(const T1 & z, T2 & y)
+{
+  assert(::pressio::ops::extent(z, 0)==::pressio::ops::extent(y, 0));
+
+  auto & y_n = impl::get_native(y);
+  const auto & z_n = impl::get_native(z);
+
+  y_n = z_n.cwiseInverse();
+}
+
+}}//end namespace pressio::ops
+#endif  // PRESSIOOPS_OPS_EIGEN_OPS_ELEMENTWISE_INVERSE_HPP_

tests/ops/ops_eigen_diag.cc

@@ -389,3 +389,26 @@ TEST(ops_eigen_diag, elementwiseMultiply)
   EXPECT_DOUBLE_EQ( y(1), 14.0);
   EXPECT_DOUBLE_EQ( y(2), 23.0);
 }
+
+TEST(ops_eigen_diag, elementwiseInverse)
+{
+  Eigen::MatrixXd M1(3,3);
+  auto y = pressio::diagonal(M1);
+
+  auto z0 = 3.;
+  auto z1 = 4.;
+  auto z2 = 5.;
+
+  Eigen::MatrixXd M2(3,3);
+  M2(0,0)=z0; M2(1,1)=z1; M2(2,2)=z2;
+  const auto z = pressio::diagonal(M2);
+
+  auto y0 = 1. / z0;
+  auto y1 = 1. / z1;
+  auto y2 = 1. / z2;
+
+  pressio::ops::elementwise_inverse(z, y);
+  EXPECT_DOUBLE_EQ( y(0), y0);
+  EXPECT_DOUBLE_EQ( y(1), y1);
+  EXPECT_DOUBLE_EQ( y(2), y2);
+}

tests/ops/ops_eigen_span.cc

@@ -266,3 +266,26 @@ TEST(ops_eigen_span, elementwiseMultiply)
   EXPECT_DOUBLE_EQ( y(1), 14.0);
   EXPECT_DOUBLE_EQ( y(2), 23.0);
 }
+
+TEST(ops_eigen_span, elementwiseInverse)
+{
+  Eigen::VectorXd M1(6);
+  auto y = pressio::span(M1,2,3);
+
+  auto z0 = 3.;
+  auto z1 = 4.;
+  auto z2 = 5.;
+
+  Eigen::VectorXd M3(6);
+  M3(2)=z0; M3(3)=z1; M3(4)=z2;
+  const auto z = pressio::span(M3,2,3);
+
+  auto y0 = 1. / z0;
+  auto y1 = 1. / z1;
+  auto y2 = 1. / z2;
+
+  pressio::ops::elementwise_inverse(z, y);
+  EXPECT_DOUBLE_EQ( y(0), y0);
+  EXPECT_DOUBLE_EQ( y(1), y1);
+  EXPECT_DOUBLE_EQ( y(2), y2);
+}

tests/ops/ops_eigen_vector.cc

@@ -414,3 +414,14 @@ TEST(ops_eigen_vector, elementwiseMultiply)
   EXPECT_DOUBLE_EQ( y(1), 12.0);
   EXPECT_DOUBLE_EQ( y(2), 20.0);
 }
+
+TEST(ops_eigen_vector, elementwiseInverse)
+{
+  V_t y; y << 0.,0.,0.;
+  V_t z; z << 3.,4.,5.;
+
+  pressio::ops::elementwise_inverse(z, y);
+  EXPECT_DOUBLE_EQ( y(0), 1./3.);
+  EXPECT_DOUBLE_EQ( y(1), 1./4.);
+  EXPECT_DOUBLE_EQ( y(2), 1./5.);
+}