diff --git a/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h b/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
index f84da9139..161e6b572 100644
--- a/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
+++ b/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
@@ -325,6 +325,22 @@ class SelfAdjointEigenSolver {
     return m_eivec * m_eivalues.cwiseSqrt().asDiagonal() * m_eivec.adjoint();
   }
 
+  /** \brief Computes the matrix exponential the matrix.
+   *
+   * \returns the matrix exponential the matrix.
+   *
+   * \pre The eigenvalues and eigenvectors of a positive-definite matrix
+   * have been computed before.
+   *
+   * \sa operatorInverseSqrt(), operatorSqrt(),
+   * <a href="unsupported/group__MatrixFunctions__Module.html">MatrixFunctions Module</a>
+   */
+  EIGEN_DEVICE_FUNC MatrixType operatorExp() const {
+    eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
+    eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
+    return m_eivec * m_eivalues.array().exp().matrix().asDiagonal() * m_eivec.adjoint();
+  }
+
   /** \brief Computes the inverse square root of the matrix.
    *
    * \returns the inverse positive-definite square root of the matrix
diff --git a/test/eigensolver_selfadjoint.cpp b/test/eigensolver_selfadjoint.cpp
index 04e727d83..12bc3b1f7 100644
--- a/test/eigensolver_selfadjoint.cpp
+++ b/test/eigensolver_selfadjoint.cpp
@@ -13,6 +13,7 @@
 #include <limits>
 #include <Eigen/Eigenvalues>
 #include <Eigen/SparseCore>
+#include <unsupported/Eigen/MatrixFunctions>
 
 template <typename MatrixType>
 void selfadjointeigensolver_essential_check(const MatrixType& m) {
@@ -135,11 +136,13 @@ void selfadjointeigensolver(const MatrixType& m) {
   VERIFY_RAISES_ASSERT(eiSymmUninitialized.eigenvectors());
   VERIFY_RAISES_ASSERT(eiSymmUninitialized.operatorSqrt());
   VERIFY_RAISES_ASSERT(eiSymmUninitialized.operatorInverseSqrt());
+  VERIFY_RAISES_ASSERT(eiSymmUninitialized.operatorExp());
 
   eiSymmUninitialized.compute(symmA, false);
   VERIFY_RAISES_ASSERT(eiSymmUninitialized.eigenvectors());
   VERIFY_RAISES_ASSERT(eiSymmUninitialized.operatorSqrt());
   VERIFY_RAISES_ASSERT(eiSymmUninitialized.operatorInverseSqrt());
+  VERIFY_RAISES_ASSERT(eiSymmUninitialized.operatorExp());
 
   // test Tridiagonalization's methods
   Tridiagonalization<MatrixType> tridiag(symmC);
@@ -167,6 +170,14 @@ void selfadjointeigensolver(const MatrixType& m) {
                                             eiSymmTridiag.eigenvectors().real().transpose());
   }
 
+  // Test matrix expponential from eigendecomposition.
+  // First scale to avoid overflow.
+  symmB = symmB / symmB.norm();
+  eiSymm.compute(symmB);
+  MatrixType expSymmB = eiSymm.operatorExp();
+  symmB = symmB.template selfadjointView<Lower>();
+  VERIFY_IS_APPROX(expSymmB, symmB.exp());
+
   if (rows > 1 && rows < 20) {
     // Test matrix with NaN
     symmC(0, 0) = std::numeric_limits<typename MatrixType::RealScalar>::quiet_NaN();