[flang] Fix crash in folding of DPROD() with non-scalar arguments

[klausler]

Summary

DPROD(x,y) is defined as DBLE(x)*DBLE(y) and that's exactly how
the implementation of its rewriting and possible folding should
be implemented, instead of the current code that only works when
both arguments are scalar and crashes otherwise.

[vzakhari]

Looks good! Thank you for the fix!

flang/lib/Evaluate/fold-real.cpp

return FoldElementalIntrinsic<T, T, T>(context, std::move(funcRef),
           if (result.flags.test(RealFlag::Overflow)) {
             context.messages().Say("DIM intrinsic folding overflow"_warn_en_US);
           }
           return result.value;
         }));
   } else if (name == "dot_product") {
     return FoldDotProduct<T>(context, std::move(funcRef));
   } else if (name == "dprod") {
-    if (auto scalars{GetScalarConstantArguments<T, T>(context, args)}) {
+    // Rewrite DPROD(x,y) -> DBLE(x)*DBLE(y)
+    if (args.at(0) && args.at(1)) {
+      const auto *xExpr{args[0]->UnwrapExpr()};
+      const auto *yExpr{args[1]->UnwrapExpr()};
+      if (xExpr && yExpr) {
         return Fold(context,
-          Expr<T>{Multiply<T>{
-              Expr<T>{std::get<0>(*scalars)}, Expr<T>{std::get<1>(*scalars)}}});
+            ToReal<T::kind>(context, common::Clone(*xExpr)) *
+                ToReal<T::kind>(context, common::Clone(*yExpr)));
+      }
     }
   } else if (name == "epsilon") {
     return Expr<T>{Scalar<T>::EPSILON()};
   } else if (name == "fraction") {
     return FoldElementalIntrinsic<T, T>(context, std::move(funcRef),
         ScalarFunc<T, T>(
             [](const Scalar<T> &x) -> Scalar<T> { return x.FRACTION(); }));
   } else if (name == "huge") {

flang/test/Evaluate/fold-dprod.f90

+! RUN: %python %S/test_folding.py %s %flang_fc1
+! Tests folding of DPROD()
+module m
+  logical, parameter :: test_kind = kind(dprod(2., 3.)) == kind(0.d0)
+  logical, parameter :: test_ss = dprod(2., 3.) == 6.d0
+  logical, parameter :: test_sv = all(dprod(2., [3.,4.]) == [6.d0,8.d0])
+  logical, parameter :: test_vv = all(dprod([2.,3.], [4.,5.]) == [8.d0,15.0d0])
+end