diff --git a/flang/include/flang/Evaluate/real.h b/flang/include/flang/Evaluate/real.h
--- a/flang/include/flang/Evaluate/real.h
+++ b/flang/include/flang/Evaluate/real.h
@@ -193,6 +193,10 @@
                 .IBSET(significandBits - 2)};
   }
 
+  static constexpr Real PositiveZero() { return Real{}; }
+
+  static constexpr Real NegativeZero() { return {Word{}.IBSET(bits - 1)}; }
+
   static constexpr Real Infinity(bool negative) {
     Word infinity{maxExponent};
     infinity = infinity.SHIFTL(significandBits);
diff --git a/flang/lib/Evaluate/real.cpp b/flang/lib/Evaluate/real.cpp
--- a/flang/lib/Evaluate/real.cpp
+++ b/flang/lib/Evaluate/real.cpp
@@ -98,7 +98,7 @@
     if (order == Ordering::Equal) {
       // x + (-x) -> +0.0 unless rounding is directed downwards
       if (rounding.mode == common::RoundingMode::Down) {
-        result.value.word_ = result.value.word_.IBSET(bits - 1); // -0.0
+        result.value = NegativeZero();
       }
       return result;
     }
@@ -221,7 +221,7 @@
       }
     } else if (IsZero() || y.IsInfinite()) { // 0/x, x/Inf -> 0
       if (isNegative) {
-        result.value.word_ = result.value.word_.IBSET(bits - 1);
+        result.value = NegativeZero();
       }
     } else {
       // dividend and divisor are both finite and nonzero numbers
@@ -272,13 +272,15 @@
   } else if (IsNegative()) {
     if (IsZero()) {
       // SQRT(-0) == -0 in IEEE-754.
-      result.value.word_ = result.value.word_.IBSET(bits - 1);
+      result.value = NegativeZero();
     } else {
       result.value = NotANumber();
     }
   } else if (IsInfinite()) {
     // SQRT(+Inf) == +Inf
     result.value = Infinity(false);
+  } else if (IsZero()) {
+    result.value = PositiveZero();
   } else {
     int expo{UnbiasedExponent()};
     if (expo < -1 || expo > 1) {
diff --git a/flang/test/Evaluate/folding28.f90 b/flang/test/Evaluate/folding28.f90
--- a/flang/test/Evaluate/folding28.f90
+++ b/flang/test/Evaluate/folding28.f90
@@ -44,4 +44,9 @@
   logical, parameter :: test_before_1 = sqrt_before_1 == before_1
   real(4), parameter :: sq_sqrt_before_1 = sqrt_before_1 * sqrt_before_1
   logical, parameter :: test_sq_before_1 = sq_sqrt_before_1 < before_1
+  ! ICE at 0.0
+  real(4), parameter :: sqrt_zero_4 = sqrt(0.0)
+  logical, parameter :: test_sqrt_zero_4 = sqrt_zero_4 == 0.0
+  real(8), parameter :: sqrt_zero_8 = sqrt(0.0)
+  logical, parameter :: test_sqrt_zero_8 = sqrt_zero_8 == 0.0
 end module