Eigen support for special matrix objects
Functions returning specialized Eigen matrices like Eigen::DiagonalMatrix and
Eigen::SelfAdjointView--which inherit from EigenBase but not
DenseBase--isn't currently allowed; such classes are explicitly copyable
into a Matrix (by definition), and so we can support functions that
return them by copying the value into a Matrix then casting that
resulting dense Matrix into a numpy.ndarray. This commit does exactly
that.
diff --git a/example/eigen.cpp b/example/eigen.cpp
index 9c33c55..c88bee2 100644
--- a/example/eigen.cpp
+++ b/example/eigen.cpp
@@ -66,6 +66,22 @@
return x.block(start_row, start_col, block_rows, block_cols);
});
+ // Returns a DiagonalMatrix with diagonal (1,2,3,...)
+ m.def("incr_diag", [](int k) {
+ Eigen::DiagonalMatrix<int, Eigen::Dynamic> m(k);
+ for (int i = 0; i < k; i++) m.diagonal()[i] = i+1;
+ return m;
+ });
+
+ // Returns a SelfAdjointView referencing the lower triangle of m
+ m.def("symmetric_lower", [](const Eigen::MatrixXi &m) {
+ return m.selfadjointView<Eigen::Lower>();
+ });
+ // Returns a SelfAdjointView referencing the lower triangle of m
+ m.def("symmetric_upper", [](const Eigen::MatrixXi &m) {
+ return m.selfadjointView<Eigen::Upper>();
+ });
+
m.def("fixed_r", [mat]() -> FixedMatrixR {
return FixedMatrixR(mat);
});
diff --git a/example/eigen.py b/example/eigen.py
index 04078b9..5f7ec51 100644
--- a/example/eigen.py
+++ b/example/eigen.py
@@ -14,6 +14,7 @@
from example import cholesky1, cholesky2, cholesky3, cholesky4, cholesky5, cholesky6
from example import diagonal, diagonal_1, diagonal_n
from example import block
+from example import incr_diag, symmetric_upper, symmetric_lower
try:
import numpy as np
import scipy
@@ -88,3 +89,20 @@
print("block(2,1,3,3) %s" % ("OK" if (block(ref, 2, 1, 3, 3) == ref[2:5, 1:4]).all() else "FAILED"))
print("block(1,4,4,2) %s" % ("OK" if (block(ref, 1, 4, 4, 2) == ref[1:, 4:]).all() else "FAILED"))
print("block(1,4,3,2) %s" % ("OK" if (block(ref, 1, 4, 3, 2) == ref[1:4, 4:]).all() else "FAILED"))
+
+print("incr_diag %s" % ("OK" if (incr_diag(7) == np.diag([1,2,3,4,5,6,7])).all() else "FAILED"))
+
+asymm = np.array([
+ [1, 2, 3, 4],
+ [5, 6, 7, 8],
+ [9, 10,11,12],
+ [13,14,15,16]])
+symm_lower = np.array(asymm)
+symm_upper = np.array(asymm)
+for i in range(4):
+ for j in range(i+1, 4):
+ symm_lower[i,j] = symm_lower[j,i]
+ symm_upper[j,i] = symm_upper[i,j]
+
+print("symmetric_lower %s" % ("OK" if (symmetric_lower(asymm) == symm_lower).all() else "FAILED"))
+print("symmetric_upper %s" % ("OK" if (symmetric_upper(asymm) == symm_upper).all() else "FAILED"))
diff --git a/example/eigen.ref b/example/eigen.ref
index 8ccd1f4..755012f 100644
--- a/example/eigen.ref
+++ b/example/eigen.ref
@@ -50,3 +50,6 @@
block(2,1,3,3) OK
block(1,4,4,2) OK
block(1,4,3,2) OK
+incr_diag OK
+symmetric_lower OK
+symmetric_upper OK