diff --git a/Eigen/src/Core/functors/UnaryFunctors.h b/Eigen/src/Core/functors/UnaryFunctors.h
index 0a6f52a72..a9ddaabb6 100644
--- a/Eigen/src/Core/functors/UnaryFunctors.h
+++ b/Eigen/src/Core/functors/UnaryFunctors.h
@@ -323,6 +323,11 @@ struct functor_traits<scalar_real_ref_op<Scalar>> {
  */
 template <typename Scalar>
 struct scalar_imag_ref_op {
+  // For a real Scalar, numext::imag_ref returns the temporary RealScalar(0) by value; binding it to the
+  // reference return types below would dangle (issue #3096). Real-valued objects get a read-only zero
+  // expression from imag() instead (see NonConstImagReturnType in CommonCwiseUnaryOps.inc).
+  static_assert(NumTraits<Scalar>::IsComplex,
+                "THE IMAGINARY PART OF A REAL-VALUED OBJECT IS NOT AN LVALUE. USE THE READ-ONLY imag() OVERLOAD.");
   typedef typename NumTraits<Scalar>::Real result_type;
   EIGEN_DEVICE_FUNC constexpr EIGEN_STRONG_INLINE result_type& operator()(Scalar& a) const {
     return numext::imag_ref(a);
diff --git a/Eigen/src/plugins/CommonCwiseUnaryOps.inc b/Eigen/src/plugins/CommonCwiseUnaryOps.inc
index 9506960e8..b2c4eb361 100644
--- a/Eigen/src/plugins/CommonCwiseUnaryOps.inc
+++ b/Eigen/src/plugins/CommonCwiseUnaryOps.inc
@@ -27,8 +27,12 @@ typedef std::conditional_t<NumTraits<Scalar>::IsComplex, CwiseUnaryView<internal
     NonConstRealReturnType;
 /** \internal the return type of imag() const */
 typedef CwiseUnaryOp<internal::scalar_imag_op<Scalar>, const Derived> ImagReturnType;
-/** \internal the return type of imag() */
-typedef CwiseUnaryView<internal::scalar_imag_ref_op<Scalar>, Derived> NonConstImagReturnType;
+/** \internal the return type of imag(). The imaginary part of a real-valued object is identically zero and not an
+ * lvalue, so for real scalars the non-const overload returns the same read-only zero expression as the const one
+ * instead of a writable view (whose functor would return a dangling reference to a temporary; see issue #3096). */
+typedef std::conditional_t<NumTraits<Scalar>::IsComplex, CwiseUnaryView<internal::scalar_imag_ref_op<Scalar>, Derived>,
+                           ImagReturnType>
+    NonConstImagReturnType;
 
 typedef CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const Derived> NegativeReturnType;
 
@@ -161,6 +165,8 @@ EIGEN_DOC_UNARY_ADDONS(real, real part function)
 EIGEN_DEVICE_FUNC constexpr inline NonConstRealReturnType real() { return NonConstRealReturnType(derived()); }
 
 /// \returns a non const expression of the imaginary part of \c *this.
+/// For real-valued objects the imaginary part is identically zero and not writable; this then returns the
+/// same read-only zero expression as the const overload.
 ///
 EIGEN_DOC_UNARY_ADDONS(imag, imaginary part function)
 ///
diff --git a/test/basicstuff.cpp b/test/basicstuff.cpp
index 94fed53a2..1ac9064d3 100644
--- a/test/basicstuff.cpp
+++ b/test/basicstuff.cpp
@@ -196,6 +196,28 @@ void basicStuffComplex(const MatrixType& m) {
   VERIFY(!static_cast<const MatrixType&>(cm).imag().isZero());
 }
 
+// Regression test for issue #3096: reading imag() of a *non-const* real-valued object used to go
+// through the writable CwiseUnaryView, whose functor returns a dangling reference to a temporary
+// for real scalars (segfault with GCC, garbage with Clang). For real scalars the non-const imag()
+// overload must return the same read-only zero expression as the const one.
+template <typename MatrixType>
+void basicStuffRealImag(const MatrixType& m) {
+  typedef typename MatrixType::RealScalar RealScalar;
+
+  Index rows = m.rows();
+  Index cols = m.cols();
+
+  MatrixType m1 = MatrixType::Random(rows, cols);
+  // Non-const lvalue.
+  VERIFY_IS_EQUAL(m1.imag().cwiseAbs().sum(), RealScalar(0));
+  // Non-const expression wrappers (temporaries), the case that used to crash.
+  VERIFY((m1.array().imag() == RealScalar(0)).all());
+  VERIFY_IS_EQUAL(m1.block(0, 0, rows, cols).imag().cwiseAbs().sum(), RealScalar(0));
+  // real() of a real-valued object aliases the object itself, both in the const and non-const case.
+  VERIFY(m1.real().data() == m1.data());
+  VERIFY(static_cast<const MatrixType&>(m1).real().data() == m1.data());
+}
+
 template <typename SrcScalar, typename TgtScalar>
 struct casting_test {
   static void run() {
@@ -347,6 +369,11 @@ EIGEN_DECLARE_TEST(basicstuff) {
         MatrixXcf(internal::random<int>(1, EIGEN_TEST_MAX_SIZE), internal::random<int>(1, EIGEN_TEST_MAX_SIZE))));
     CALL_SUBTEST_5(basicStuffComplex(
         MatrixXcd(internal::random<int>(1, EIGEN_TEST_MAX_SIZE), internal::random<int>(1, EIGEN_TEST_MAX_SIZE))));
+
+    CALL_SUBTEST_1(basicStuffRealImag(Matrix<float, 1, 1>()));
+    CALL_SUBTEST_2(basicStuffRealImag(Matrix4d()));
+    CALL_SUBTEST_7(basicStuffRealImag(
+        MatrixXd(internal::random<int>(1, EIGEN_TEST_MAX_SIZE), internal::random<int>(1, EIGEN_TEST_MAX_SIZE))));
   }
 
   CALL_SUBTEST_1(fixedSizeMatrixConstruction<unsigned char>());