1 //===-- runtime/dot-product.cpp -------------------------------------------===//
2 //
3 // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
4 // See https://llvm.org/LICENSE.txt for license information.
5 // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
6 //
7 //===----------------------------------------------------------------------===//
8 
9 #include "cpp-type.h"
10 #include "descriptor.h"
11 #include "reduction.h"
12 #include "terminator.h"
13 #include "tools.h"
14 #include <cinttypes>
15 
16 namespace Fortran::runtime {
17 
18 template <typename ACCUMULATOR>
19 static inline auto DoDotProduct(const Descriptor &x, const Descriptor &y,
20     Terminator &terminator) -> typename ACCUMULATOR::Result {
21   RUNTIME_CHECK(terminator, x.rank() == 1 && y.rank() == 1);
22   SubscriptValue n{x.GetDimension(0).Extent()};
23   if (SubscriptValue yN{y.GetDimension(0).Extent()}; yN != n) {
24     terminator.Crash(
25         "DOT_PRODUCT: SIZE(VECTOR_A) is %jd but SIZE(VECTOR_B) is %jd",
26         static_cast<std::intmax_t>(n), static_cast<std::intmax_t>(yN));
27   }
28   SubscriptValue xAt{x.GetDimension(0).LowerBound()};
29   SubscriptValue yAt{y.GetDimension(0).LowerBound()};
30   ACCUMULATOR accumulator{x, y};
31   for (SubscriptValue j{0}; j < n; ++j) {
32     accumulator.Accumulate(xAt++, yAt++);
33   }
34   return accumulator.GetResult();
35 }
36 
37 template <TypeCategory RCAT, int RKIND,
38     template <typename, TypeCategory, typename, typename> class ACCUM>
39 struct DotProduct {
40   using Result = CppTypeFor<RCAT, RKIND>;
41   template <TypeCategory XCAT, int XKIND> struct DP1 {
42     template <TypeCategory YCAT, int YKIND> struct DP2 {
43       Result operator()(const Descriptor &x, const Descriptor &y,
44           Terminator &terminator) const {
45         if constexpr (constexpr auto resultType{
46                           GetResultType(XCAT, XKIND, YCAT, YKIND)}) {
47           if constexpr (resultType->first == RCAT &&
48               resultType->second <= RKIND) {
49             using Accum = ACCUM<Result, XCAT, CppTypeFor<XCAT, XKIND>,
50                 CppTypeFor<YCAT, YKIND>>;
51             return DoDotProduct<Accum>(x, y, terminator);
52           }
53         }
54         terminator.Crash(
55             "DOT_PRODUCT(%d(%d)): bad operand types (%d(%d), %d(%d))",
56             static_cast<int>(RCAT), RKIND, static_cast<int>(XCAT), XKIND,
57             static_cast<int>(YCAT), YKIND);
58       }
59     };
60     Result operator()(const Descriptor &x, const Descriptor &y,
61         Terminator &terminator, TypeCategory yCat, int yKind) const {
62       return ApplyType<DP2, Result>(yCat, yKind, terminator, x, y, terminator);
63     }
64   };
65   Result operator()(const Descriptor &x, const Descriptor &y,
66       const char *source, int line) const {
67     Terminator terminator{source, line};
68     auto xCatKind{x.type().GetCategoryAndKind()};
69     auto yCatKind{y.type().GetCategoryAndKind()};
70     RUNTIME_CHECK(terminator, xCatKind.has_value() && yCatKind.has_value());
71     return ApplyType<DP1, Result>(xCatKind->first, xCatKind->second, terminator,
72         x, y, terminator, yCatKind->first, yCatKind->second);
73   }
74 };
75 
76 template <typename RESULT, TypeCategory XCAT, typename XT, typename YT>
77 class NumericAccumulator {
78 public:
79   using Result = RESULT;
80   NumericAccumulator(const Descriptor &x, const Descriptor &y) : x_{x}, y_{y} {}
81   void Accumulate(SubscriptValue xAt, SubscriptValue yAt) {
82     if constexpr (XCAT == TypeCategory::Complex) {
83       sum_ += std::conj(static_cast<Result>(*x_.Element<XT>(&xAt))) *
84           static_cast<Result>(*y_.Element<YT>(&yAt));
85     } else {
86       sum_ += static_cast<Result>(*x_.Element<XT>(&xAt)) *
87           static_cast<Result>(*y_.Element<YT>(&yAt));
88     }
89   }
90   Result GetResult() const { return sum_; }
91 
92 private:
93   const Descriptor &x_, &y_;
94   Result sum_{0};
95 };
96 
97 template <typename, TypeCategory, typename XT, typename YT>
98 class LogicalAccumulator {
99 public:
100   using Result = bool;
101   LogicalAccumulator(const Descriptor &x, const Descriptor &y) : x_{x}, y_{y} {}
102   void Accumulate(SubscriptValue xAt, SubscriptValue yAt) {
103     result_ = result_ ||
104         (IsLogicalElementTrue(x_, &xAt) && IsLogicalElementTrue(y_, &yAt));
105   }
106   bool GetResult() const { return result_; }
107 
108 private:
109   const Descriptor &x_, &y_;
110   bool result_{false};
111 };
112 
113 extern "C" {
114 std::int8_t RTNAME(DotProductInteger1)(
115     const Descriptor &x, const Descriptor &y, const char *source, int line) {
116   return DotProduct<TypeCategory::Integer, 8, NumericAccumulator>{}(
117       x, y, source, line);
118 }
119 std::int16_t RTNAME(DotProductInteger2)(
120     const Descriptor &x, const Descriptor &y, const char *source, int line) {
121   return DotProduct<TypeCategory::Integer, 8, NumericAccumulator>{}(
122       x, y, source, line);
123 }
124 std::int32_t RTNAME(DotProductInteger4)(
125     const Descriptor &x, const Descriptor &y, const char *source, int line) {
126   return DotProduct<TypeCategory::Integer, 8, NumericAccumulator>{}(
127       x, y, source, line);
128 }
129 std::int64_t RTNAME(DotProductInteger8)(
130     const Descriptor &x, const Descriptor &y, const char *source, int line) {
131   return DotProduct<TypeCategory::Integer, 8, NumericAccumulator>{}(
132       x, y, source, line);
133 }
134 #ifdef __SIZEOF_INT128__
135 common::int128_t RTNAME(DotProductInteger16)(
136     const Descriptor &x, const Descriptor &y, const char *source, int line) {
137   return DotProduct<TypeCategory::Integer, 16, NumericAccumulator>{}(
138       x, y, source, line);
139 }
140 #endif
141 
142 // TODO: REAL/COMPLEX(2 & 3)
143 float RTNAME(DotProductReal4)(
144     const Descriptor &x, const Descriptor &y, const char *source, int line) {
145   return DotProduct<TypeCategory::Real, 8, NumericAccumulator>{}(
146       x, y, source, line);
147 }
148 double RTNAME(DotProductReal8)(
149     const Descriptor &x, const Descriptor &y, const char *source, int line) {
150   return DotProduct<TypeCategory::Real, 8, NumericAccumulator>{}(
151       x, y, source, line);
152 }
153 #if LONG_DOUBLE == 80
154 long double RTNAME(DotProductReal10)(
155     const Descriptor &x, const Descriptor &y, const char *source, int line) {
156   return DotProduct<TypeCategory::Real, 10, NumericAccumulator>{}(
157       x, y, source, line);
158 }
159 #elif LONG_DOUBLE == 128
160 long double RTNAME(DotProductReal16)(
161     const Descriptor &x, const Descriptor &y, const char *source, int line) {
162   return DotProduct<TypeCategory::Real, 16, NumericAccumulator>{}(
163       x, y, source, line);
164 }
165 #endif
166 
167 void RTNAME(CppDotProductComplex4)(std::complex<float> &result,
168     const Descriptor &x, const Descriptor &y, const char *source, int line) {
169   auto z{DotProduct<TypeCategory::Complex, 8, NumericAccumulator>{}(
170       x, y, source, line)};
171   result = std::complex<float>{
172       static_cast<float>(z.real()), static_cast<float>(z.imag())};
173 }
174 void RTNAME(CppDotProductComplex8)(std::complex<double> &result,
175     const Descriptor &x, const Descriptor &y, const char *source, int line) {
176   result = DotProduct<TypeCategory::Complex, 8, NumericAccumulator>{}(
177       x, y, source, line);
178 }
179 #if LONG_DOUBLE == 80
180 void RTNAME(CppDotProductComplex10)(std::complex<long double> &result,
181     const Descriptor &x, const Descriptor &y, const char *source, int line) {
182   result = DotProduct<TypeCategory::Complex, 10, NumericAccumulator>{}(
183       x, y, source, line);
184 }
185 #elif LONG_DOUBLE == 128
186 void RTNAME(CppDotProductComplex16)(std::complex<long double> &result,
187     const Descriptor &x, const Descriptor &y, const char *source, int line) {
188   result = DotProduct<TypeCategory::Complex, 16, NumericAccumulator>{}(
189       x, y, source, line);
190 }
191 #endif
192 
193 bool RTNAME(DotProductLogical)(
194     const Descriptor &x, const Descriptor &y, const char *source, int line) {
195   return DotProduct<TypeCategory::Logical, 1, LogicalAccumulator>{}(
196       x, y, source, line);
197 }
198 } // extern "C"
199 } // namespace Fortran::runtime
200