1 //===-- runtime/extrema.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 // Implements MAXLOC, MINLOC, MAXVAL, & MINVAL for all required operand types
10 // and shapes and (for MAXLOC & MINLOC) result integer kinds.  Also implements
11 // NORM2 using common infrastructure.
12 
13 #include "reduction-templates.h"
14 #include "flang/Common/long-double.h"
15 #include "flang/Runtime/character.h"
16 #include "flang/Runtime/reduction.h"
17 #include <algorithm>
18 #include <cinttypes>
19 #include <cmath>
20 #include <optional>
21 
22 namespace Fortran::runtime {
23 
24 // MAXLOC & MINLOC
25 
26 template <typename T, bool IS_MAX, bool BACK> struct NumericCompare {
27   using Type = T;
28   explicit NumericCompare(std::size_t /*elemLen; ignored*/) {}
29   bool operator()(const T &value, const T &previous) const {
30     if (value == previous) {
31       return BACK;
32     } else if constexpr (IS_MAX) {
33       return value > previous;
34     } else {
35       return value < previous;
36     }
37   }
38 };
39 
40 template <typename T, bool IS_MAX, bool BACK> class CharacterCompare {
41 public:
42   using Type = T;
43   explicit CharacterCompare(std::size_t elemLen)
44       : chars_{elemLen / sizeof(T)} {}
45   bool operator()(const T &value, const T &previous) const {
46     int cmp{CharacterScalarCompare<T>(&value, &previous, chars_, chars_)};
47     if (cmp == 0) {
48       return BACK;
49     } else if constexpr (IS_MAX) {
50       return cmp > 0;
51     } else {
52       return cmp < 0;
53     }
54   }
55 
56 private:
57   std::size_t chars_;
58 };
59 
60 template <typename COMPARE> class ExtremumLocAccumulator {
61 public:
62   using Type = typename COMPARE::Type;
63   ExtremumLocAccumulator(const Descriptor &array, std::size_t chars = 0)
64       : array_{array}, argRank_{array.rank()}, compare_{array.ElementBytes()} {
65     Reinitialize();
66   }
67   void Reinitialize() {
68     // per standard: result indices are all zero if no data
69     for (int j{0}; j < argRank_; ++j) {
70       extremumLoc_[j] = 0;
71     }
72     previous_ = nullptr;
73   }
74   int argRank() const { return argRank_; }
75   template <typename A> void GetResult(A *p, int zeroBasedDim = -1) {
76     if (zeroBasedDim >= 0) {
77       *p = extremumLoc_[zeroBasedDim] -
78           array_.GetDimension(zeroBasedDim).LowerBound() + 1;
79     } else {
80       for (int j{0}; j < argRank_; ++j) {
81         p[j] = extremumLoc_[j] - array_.GetDimension(j).LowerBound() + 1;
82       }
83     }
84   }
85   template <typename IGNORED> bool AccumulateAt(const SubscriptValue at[]) {
86     const auto &value{*array_.Element<Type>(at)};
87     if (!previous_ || compare_(value, *previous_)) {
88       previous_ = &value;
89       for (int j{0}; j < argRank_; ++j) {
90         extremumLoc_[j] = at[j];
91       }
92     }
93     return true;
94   }
95 
96 private:
97   const Descriptor &array_;
98   int argRank_;
99   SubscriptValue extremumLoc_[maxRank];
100   const Type *previous_{nullptr};
101   COMPARE compare_;
102 };
103 
104 template <typename ACCUMULATOR, typename CPPTYPE>
105 static void LocationHelper(const char *intrinsic, Descriptor &result,
106     const Descriptor &x, int kind, const Descriptor *mask,
107     Terminator &terminator) {
108   ACCUMULATOR accumulator{x};
109   DoTotalReduction<CPPTYPE>(x, 0, mask, accumulator, intrinsic, terminator);
110   ApplyIntegerKind<LocationResultHelper<ACCUMULATOR>::template Functor, void>(
111       kind, terminator, accumulator, result);
112 }
113 
114 template <TypeCategory CAT, int KIND, bool IS_MAX,
115     template <typename, bool, bool> class COMPARE>
116 inline void DoMaxOrMinLoc(const char *intrinsic, Descriptor &result,
117     const Descriptor &x, int kind, const char *source, int line,
118     const Descriptor *mask, bool back) {
119   using CppType = CppTypeFor<CAT, KIND>;
120   Terminator terminator{source, line};
121   if (back) {
122     LocationHelper<ExtremumLocAccumulator<COMPARE<CppType, IS_MAX, true>>,
123         CppType>(intrinsic, result, x, kind, mask, terminator);
124   } else {
125     LocationHelper<ExtremumLocAccumulator<COMPARE<CppType, IS_MAX, false>>,
126         CppType>(intrinsic, result, x, kind, mask, terminator);
127   }
128 }
129 
130 template <TypeCategory CAT, bool IS_MAX> struct TypedMaxOrMinLocHelper {
131   template <int KIND> struct Functor {
132     void operator()(const char *intrinsic, Descriptor &result,
133         const Descriptor &x, int kind, const char *source, int line,
134         const Descriptor *mask, bool back) const {
135       DoMaxOrMinLoc<TypeCategory::Integer, KIND, IS_MAX, NumericCompare>(
136           intrinsic, result, x, kind, source, line, mask, back);
137     }
138   };
139 };
140 
141 template <bool IS_MAX>
142 inline void TypedMaxOrMinLoc(const char *intrinsic, Descriptor &result,
143     const Descriptor &x, int kind, const char *source, int line,
144     const Descriptor *mask, bool back) {
145   int rank{x.rank()};
146   SubscriptValue extent[1]{rank};
147   result.Establish(TypeCategory::Integer, kind, nullptr, 1, extent,
148       CFI_attribute_allocatable);
149   result.GetDimension(0).SetBounds(1, extent[0]);
150   Terminator terminator{source, line};
151   if (int stat{result.Allocate()}) {
152     terminator.Crash(
153         "%s: could not allocate memory for result; STAT=%d", intrinsic, stat);
154   }
155   CheckIntegerKind(terminator, kind, intrinsic);
156   auto catKind{x.type().GetCategoryAndKind()};
157   RUNTIME_CHECK(terminator, catKind.has_value());
158   switch (catKind->first) {
159   case TypeCategory::Integer:
160     ApplyIntegerKind<
161         TypedMaxOrMinLocHelper<TypeCategory::Integer, IS_MAX>::template Functor,
162         void>(catKind->second, terminator, intrinsic, result, x, kind, source,
163         line, mask, back);
164     break;
165   case TypeCategory::Real:
166     ApplyFloatingPointKind<
167         TypedMaxOrMinLocHelper<TypeCategory::Real, IS_MAX>::template Functor,
168         void>(catKind->second, terminator, intrinsic, result, x, kind, source,
169         line, mask, back);
170     break;
171   case TypeCategory::Character:
172     ApplyCharacterKind<TypedMaxOrMinLocHelper<TypeCategory::Character,
173                            IS_MAX>::template Functor,
174         void>(catKind->second, terminator, intrinsic, result, x, kind, source,
175         line, mask, back);
176     break;
177   default:
178     terminator.Crash(
179         "%s: Bad data type code (%d) for array", intrinsic, x.type().raw());
180   }
181 }
182 
183 extern "C" {
184 void RTNAME(Maxloc)(Descriptor &result, const Descriptor &x, int kind,
185     const char *source, int line, const Descriptor *mask, bool back) {
186   TypedMaxOrMinLoc<true>("MAXLOC", result, x, kind, source, line, mask, back);
187 }
188 void RTNAME(Minloc)(Descriptor &result, const Descriptor &x, int kind,
189     const char *source, int line, const Descriptor *mask, bool back) {
190   TypedMaxOrMinLoc<false>("MINLOC", result, x, kind, source, line, mask, back);
191 }
192 } // extern "C"
193 
194 // MAXLOC/MINLOC with DIM=
195 
196 template <TypeCategory CAT, int KIND, bool IS_MAX,
197     template <typename, bool, bool> class COMPARE, bool BACK>
198 static void DoPartialMaxOrMinLocDirection(const char *intrinsic,
199     Descriptor &result, const Descriptor &x, int kind, int dim,
200     const Descriptor *mask, Terminator &terminator) {
201   using CppType = CppTypeFor<CAT, KIND>;
202   using Accumulator = ExtremumLocAccumulator<COMPARE<CppType, IS_MAX, BACK>>;
203   Accumulator accumulator{x};
204   ApplyIntegerKind<PartialLocationHelper<Accumulator>::template Functor, void>(
205       kind, terminator, result, x, dim, mask, terminator, intrinsic,
206       accumulator);
207 }
208 
209 template <TypeCategory CAT, int KIND, bool IS_MAX,
210     template <typename, bool, bool> class COMPARE>
211 inline void DoPartialMaxOrMinLoc(const char *intrinsic, Descriptor &result,
212     const Descriptor &x, int kind, int dim, const Descriptor *mask, bool back,
213     Terminator &terminator) {
214   if (back) {
215     DoPartialMaxOrMinLocDirection<CAT, KIND, IS_MAX, COMPARE, true>(
216         intrinsic, result, x, kind, dim, mask, terminator);
217   } else {
218     DoPartialMaxOrMinLocDirection<CAT, KIND, IS_MAX, COMPARE, false>(
219         intrinsic, result, x, kind, dim, mask, terminator);
220   }
221 }
222 
223 template <TypeCategory CAT, bool IS_MAX,
224     template <typename, bool, bool> class COMPARE>
225 struct DoPartialMaxOrMinLocHelper {
226   template <int KIND> struct Functor {
227     void operator()(const char *intrinsic, Descriptor &result,
228         const Descriptor &x, int kind, int dim, const Descriptor *mask,
229         bool back, Terminator &terminator) const {
230       DoPartialMaxOrMinLoc<CAT, KIND, IS_MAX, COMPARE>(
231           intrinsic, result, x, kind, dim, mask, back, terminator);
232     }
233   };
234 };
235 
236 template <bool IS_MAX>
237 inline void TypedPartialMaxOrMinLoc(const char *intrinsic, Descriptor &result,
238     const Descriptor &x, int kind, int dim, const char *source, int line,
239     const Descriptor *mask, bool back) {
240   Terminator terminator{source, line};
241   CheckIntegerKind(terminator, kind, intrinsic);
242   auto catKind{x.type().GetCategoryAndKind()};
243   RUNTIME_CHECK(terminator, catKind.has_value());
244   switch (catKind->first) {
245   case TypeCategory::Integer:
246     ApplyIntegerKind<DoPartialMaxOrMinLocHelper<TypeCategory::Integer, IS_MAX,
247                          NumericCompare>::template Functor,
248         void>(catKind->second, terminator, intrinsic, result, x, kind, dim,
249         mask, back, terminator);
250     break;
251   case TypeCategory::Real:
252     ApplyFloatingPointKind<DoPartialMaxOrMinLocHelper<TypeCategory::Real,
253                                IS_MAX, NumericCompare>::template Functor,
254         void>(catKind->second, terminator, intrinsic, result, x, kind, dim,
255         mask, back, terminator);
256     break;
257   case TypeCategory::Character:
258     ApplyCharacterKind<DoPartialMaxOrMinLocHelper<TypeCategory::Character,
259                            IS_MAX, CharacterCompare>::template Functor,
260         void>(catKind->second, terminator, intrinsic, result, x, kind, dim,
261         mask, back, terminator);
262     break;
263   default:
264     terminator.Crash(
265         "%s: Bad data type code (%d) for array", intrinsic, x.type().raw());
266   }
267 }
268 
269 extern "C" {
270 void RTNAME(MaxlocDim)(Descriptor &result, const Descriptor &x, int kind,
271     int dim, const char *source, int line, const Descriptor *mask, bool back) {
272   TypedPartialMaxOrMinLoc<true>(
273       "MAXLOC", result, x, kind, dim, source, line, mask, back);
274 }
275 void RTNAME(MinlocDim)(Descriptor &result, const Descriptor &x, int kind,
276     int dim, const char *source, int line, const Descriptor *mask, bool back) {
277   TypedPartialMaxOrMinLoc<false>(
278       "MINLOC", result, x, kind, dim, source, line, mask, back);
279 }
280 } // extern "C"
281 
282 // MAXVAL and MINVAL
283 
284 template <TypeCategory CAT, int KIND, bool IS_MAXVAL> struct MaxOrMinIdentity {
285   using Type = CppTypeFor<CAT, KIND>;
286   static constexpr Type Value() {
287     return IS_MAXVAL ? std::numeric_limits<Type>::lowest()
288                      : std::numeric_limits<Type>::max();
289   }
290 };
291 
292 // std::numeric_limits<> may not know int128_t
293 template <bool IS_MAXVAL>
294 struct MaxOrMinIdentity<TypeCategory::Integer, 16, IS_MAXVAL> {
295   using Type = CppTypeFor<TypeCategory::Integer, 16>;
296   static constexpr Type Value() {
297     return IS_MAXVAL ? Type{1} << 127 : ~Type{0} >> 1;
298   }
299 };
300 
301 template <TypeCategory CAT, int KIND, bool IS_MAXVAL>
302 class NumericExtremumAccumulator {
303 public:
304   using Type = CppTypeFor<CAT, KIND>;
305   explicit NumericExtremumAccumulator(const Descriptor &array)
306       : array_{array} {}
307   void Reinitialize() {
308     extremum_ = MaxOrMinIdentity<CAT, KIND, IS_MAXVAL>::Value();
309   }
310   template <typename A> void GetResult(A *p, int /*zeroBasedDim*/ = -1) const {
311     *p = extremum_;
312   }
313   bool Accumulate(Type x) {
314     if constexpr (IS_MAXVAL) {
315       if (x > extremum_) {
316         extremum_ = x;
317       }
318     } else if (x < extremum_) {
319       extremum_ = x;
320     }
321     return true;
322   }
323   template <typename A> bool AccumulateAt(const SubscriptValue at[]) {
324     return Accumulate(*array_.Element<A>(at));
325   }
326 
327 private:
328   const Descriptor &array_;
329   Type extremum_{MaxOrMinIdentity<CAT, KIND, IS_MAXVAL>::Value()};
330 };
331 
332 template <TypeCategory CAT, int KIND, bool IS_MAXVAL>
333 inline CppTypeFor<CAT, KIND> TotalNumericMaxOrMin(const Descriptor &x,
334     const char *source, int line, int dim, const Descriptor *mask,
335     const char *intrinsic) {
336   return GetTotalReduction<CAT, KIND>(x, source, line, dim, mask,
337       NumericExtremumAccumulator<CAT, KIND, IS_MAXVAL>{x}, intrinsic);
338 }
339 
340 template <TypeCategory CAT, int KIND, typename ACCUMULATOR>
341 static void DoMaxMinNorm2(Descriptor &result, const Descriptor &x, int dim,
342     const Descriptor *mask, const char *intrinsic, Terminator &terminator) {
343   using Type = CppTypeFor<CAT, KIND>;
344   ACCUMULATOR accumulator{x};
345   if (dim == 0 || x.rank() == 1) {
346     // Total reduction
347     result.Establish(x.type(), x.ElementBytes(), nullptr, 0, nullptr,
348         CFI_attribute_allocatable);
349     if (int stat{result.Allocate()}) {
350       terminator.Crash(
351           "%s: could not allocate memory for result; STAT=%d", intrinsic, stat);
352     }
353     DoTotalReduction<Type>(x, dim, mask, accumulator, intrinsic, terminator);
354     accumulator.GetResult(result.OffsetElement<Type>());
355   } else {
356     // Partial reduction
357     PartialReduction<ACCUMULATOR, CAT, KIND>(
358         result, x, dim, mask, terminator, intrinsic, accumulator);
359   }
360 }
361 
362 template <TypeCategory CAT, bool IS_MAXVAL> struct MaxOrMinHelper {
363   template <int KIND> struct Functor {
364     void operator()(Descriptor &result, const Descriptor &x, int dim,
365         const Descriptor *mask, const char *intrinsic,
366         Terminator &terminator) const {
367       DoMaxMinNorm2<CAT, KIND,
368           NumericExtremumAccumulator<CAT, KIND, IS_MAXVAL>>(
369           result, x, dim, mask, intrinsic, terminator);
370     }
371   };
372 };
373 
374 template <bool IS_MAXVAL>
375 inline void NumericMaxOrMin(Descriptor &result, const Descriptor &x, int dim,
376     const char *source, int line, const Descriptor *mask,
377     const char *intrinsic) {
378   Terminator terminator{source, line};
379   auto type{x.type().GetCategoryAndKind()};
380   RUNTIME_CHECK(terminator, type);
381   switch (type->first) {
382   case TypeCategory::Integer:
383     ApplyIntegerKind<
384         MaxOrMinHelper<TypeCategory::Integer, IS_MAXVAL>::template Functor,
385         void>(
386         type->second, terminator, result, x, dim, mask, intrinsic, terminator);
387     break;
388   case TypeCategory::Real:
389     ApplyFloatingPointKind<
390         MaxOrMinHelper<TypeCategory::Real, IS_MAXVAL>::template Functor, void>(
391         type->second, terminator, result, x, dim, mask, intrinsic, terminator);
392     break;
393   default:
394     terminator.Crash("%s: bad type code %d", intrinsic, x.type().raw());
395   }
396 }
397 
398 template <int KIND, bool IS_MAXVAL> class CharacterExtremumAccumulator {
399 public:
400   using Type = CppTypeFor<TypeCategory::Character, KIND>;
401   explicit CharacterExtremumAccumulator(const Descriptor &array)
402       : array_{array}, charLen_{array_.ElementBytes() / KIND} {}
403   void Reinitialize() { extremum_ = nullptr; }
404   template <typename A> void GetResult(A *p, int /*zeroBasedDim*/ = -1) const {
405     static_assert(std::is_same_v<A, Type>);
406     if (extremum_) {
407       std::memcpy(p, extremum_, charLen_);
408     } else {
409       // empty array: result is all zero-valued characters
410       std::memset(p, 0, charLen_);
411     }
412   }
413   bool Accumulate(const Type *x) {
414     if (!extremum_) {
415       extremum_ = x;
416     } else {
417       int cmp{CharacterScalarCompare(x, extremum_, charLen_, charLen_)};
418       if (IS_MAXVAL == (cmp > 0)) {
419         extremum_ = x;
420       }
421     }
422     return true;
423   }
424   template <typename A> bool AccumulateAt(const SubscriptValue at[]) {
425     return Accumulate(array_.Element<A>(at));
426   }
427 
428 private:
429   const Descriptor &array_;
430   std::size_t charLen_;
431   const Type *extremum_{nullptr};
432 };
433 
434 template <bool IS_MAXVAL> struct CharacterMaxOrMinHelper {
435   template <int KIND> struct Functor {
436     void operator()(Descriptor &result, const Descriptor &x, int dim,
437         const Descriptor *mask, const char *intrinsic,
438         Terminator &terminator) const {
439       DoMaxMinNorm2<TypeCategory::Character, KIND,
440           CharacterExtremumAccumulator<KIND, IS_MAXVAL>>(
441           result, x, dim, mask, intrinsic, terminator);
442     }
443   };
444 };
445 
446 template <bool IS_MAXVAL>
447 inline void CharacterMaxOrMin(Descriptor &result, const Descriptor &x, int dim,
448     const char *source, int line, const Descriptor *mask,
449     const char *intrinsic) {
450   Terminator terminator{source, line};
451   auto type{x.type().GetCategoryAndKind()};
452   RUNTIME_CHECK(terminator, type && type->first == TypeCategory::Character);
453   ApplyCharacterKind<CharacterMaxOrMinHelper<IS_MAXVAL>::template Functor,
454       void>(
455       type->second, terminator, result, x, dim, mask, intrinsic, terminator);
456 }
457 
458 extern "C" {
459 CppTypeFor<TypeCategory::Integer, 1> RTNAME(MaxvalInteger1)(const Descriptor &x,
460     const char *source, int line, int dim, const Descriptor *mask) {
461   return TotalNumericMaxOrMin<TypeCategory::Integer, 1, true>(
462       x, source, line, dim, mask, "MAXVAL");
463 }
464 CppTypeFor<TypeCategory::Integer, 2> RTNAME(MaxvalInteger2)(const Descriptor &x,
465     const char *source, int line, int dim, const Descriptor *mask) {
466   return TotalNumericMaxOrMin<TypeCategory::Integer, 2, true>(
467       x, source, line, dim, mask, "MAXVAL");
468 }
469 CppTypeFor<TypeCategory::Integer, 4> RTNAME(MaxvalInteger4)(const Descriptor &x,
470     const char *source, int line, int dim, const Descriptor *mask) {
471   return TotalNumericMaxOrMin<TypeCategory::Integer, 4, true>(
472       x, source, line, dim, mask, "MAXVAL");
473 }
474 CppTypeFor<TypeCategory::Integer, 8> RTNAME(MaxvalInteger8)(const Descriptor &x,
475     const char *source, int line, int dim, const Descriptor *mask) {
476   return TotalNumericMaxOrMin<TypeCategory::Integer, 8, true>(
477       x, source, line, dim, mask, "MAXVAL");
478 }
479 #ifdef __SIZEOF_INT128__
480 CppTypeFor<TypeCategory::Integer, 16> RTNAME(MaxvalInteger16)(
481     const Descriptor &x, const char *source, int line, int dim,
482     const Descriptor *mask) {
483   return TotalNumericMaxOrMin<TypeCategory::Integer, 16, true>(
484       x, source, line, dim, mask, "MAXVAL");
485 }
486 #endif
487 
488 // TODO: REAL(2 & 3)
489 CppTypeFor<TypeCategory::Real, 4> RTNAME(MaxvalReal4)(const Descriptor &x,
490     const char *source, int line, int dim, const Descriptor *mask) {
491   return TotalNumericMaxOrMin<TypeCategory::Real, 4, true>(
492       x, source, line, dim, mask, "MAXVAL");
493 }
494 CppTypeFor<TypeCategory::Real, 8> RTNAME(MaxvalReal8)(const Descriptor &x,
495     const char *source, int line, int dim, const Descriptor *mask) {
496   return TotalNumericMaxOrMin<TypeCategory::Real, 8, true>(
497       x, source, line, dim, mask, "MAXVAL");
498 }
499 #if LONG_DOUBLE == 80
500 CppTypeFor<TypeCategory::Real, 10> RTNAME(MaxvalReal10)(const Descriptor &x,
501     const char *source, int line, int dim, const Descriptor *mask) {
502   return TotalNumericMaxOrMin<TypeCategory::Real, 10, true>(
503       x, source, line, dim, mask, "MAXVAL");
504 }
505 #elif LONG_DOUBLE == 128
506 CppTypeFor<TypeCategory::Real, 16> RTNAME(MaxvalReal16)(const Descriptor &x,
507     const char *source, int line, int dim, const Descriptor *mask) {
508   return TotalNumericMaxOrMin<TypeCategory::Real, 16, true>(
509       x, source, line, dim, mask, "MAXVAL");
510 }
511 #endif
512 
513 void RTNAME(MaxvalCharacter)(Descriptor &result, const Descriptor &x,
514     const char *source, int line, const Descriptor *mask) {
515   CharacterMaxOrMin<true>(result, x, 0, source, line, mask, "MAXVAL");
516 }
517 
518 CppTypeFor<TypeCategory::Integer, 1> RTNAME(MinvalInteger1)(const Descriptor &x,
519     const char *source, int line, int dim, const Descriptor *mask) {
520   return TotalNumericMaxOrMin<TypeCategory::Integer, 1, false>(
521       x, source, line, dim, mask, "MINVAL");
522 }
523 CppTypeFor<TypeCategory::Integer, 2> RTNAME(MinvalInteger2)(const Descriptor &x,
524     const char *source, int line, int dim, const Descriptor *mask) {
525   return TotalNumericMaxOrMin<TypeCategory::Integer, 2, false>(
526       x, source, line, dim, mask, "MINVAL");
527 }
528 CppTypeFor<TypeCategory::Integer, 4> RTNAME(MinvalInteger4)(const Descriptor &x,
529     const char *source, int line, int dim, const Descriptor *mask) {
530   return TotalNumericMaxOrMin<TypeCategory::Integer, 4, false>(
531       x, source, line, dim, mask, "MINVAL");
532 }
533 CppTypeFor<TypeCategory::Integer, 8> RTNAME(MinvalInteger8)(const Descriptor &x,
534     const char *source, int line, int dim, const Descriptor *mask) {
535   return TotalNumericMaxOrMin<TypeCategory::Integer, 8, false>(
536       x, source, line, dim, mask, "MINVAL");
537 }
538 #ifdef __SIZEOF_INT128__
539 CppTypeFor<TypeCategory::Integer, 16> RTNAME(MinvalInteger16)(
540     const Descriptor &x, const char *source, int line, int dim,
541     const Descriptor *mask) {
542   return TotalNumericMaxOrMin<TypeCategory::Integer, 16, false>(
543       x, source, line, dim, mask, "MINVAL");
544 }
545 #endif
546 
547 // TODO: REAL(2 & 3)
548 CppTypeFor<TypeCategory::Real, 4> RTNAME(MinvalReal4)(const Descriptor &x,
549     const char *source, int line, int dim, const Descriptor *mask) {
550   return TotalNumericMaxOrMin<TypeCategory::Real, 4, false>(
551       x, source, line, dim, mask, "MINVAL");
552 }
553 CppTypeFor<TypeCategory::Real, 8> RTNAME(MinvalReal8)(const Descriptor &x,
554     const char *source, int line, int dim, const Descriptor *mask) {
555   return TotalNumericMaxOrMin<TypeCategory::Real, 8, false>(
556       x, source, line, dim, mask, "MINVAL");
557 }
558 #if LONG_DOUBLE == 80
559 CppTypeFor<TypeCategory::Real, 10> RTNAME(MinvalReal10)(const Descriptor &x,
560     const char *source, int line, int dim, const Descriptor *mask) {
561   return TotalNumericMaxOrMin<TypeCategory::Real, 10, false>(
562       x, source, line, dim, mask, "MINVAL");
563 }
564 #elif LONG_DOUBLE == 128
565 CppTypeFor<TypeCategory::Real, 16> RTNAME(MinvalReal16)(const Descriptor &x,
566     const char *source, int line, int dim, const Descriptor *mask) {
567   return TotalNumericMaxOrMin<TypeCategory::Real, 16, false>(
568       x, source, line, dim, mask, "MINVAL");
569 }
570 #endif
571 
572 void RTNAME(MinvalCharacter)(Descriptor &result, const Descriptor &x,
573     const char *source, int line, const Descriptor *mask) {
574   CharacterMaxOrMin<false>(result, x, 0, source, line, mask, "MINVAL");
575 }
576 
577 void RTNAME(MaxvalDim)(Descriptor &result, const Descriptor &x, int dim,
578     const char *source, int line, const Descriptor *mask) {
579   if (x.type().IsCharacter()) {
580     CharacterMaxOrMin<true>(result, x, dim, source, line, mask, "MAXVAL");
581   } else {
582     NumericMaxOrMin<true>(result, x, dim, source, line, mask, "MAXVAL");
583   }
584 }
585 void RTNAME(MinvalDim)(Descriptor &result, const Descriptor &x, int dim,
586     const char *source, int line, const Descriptor *mask) {
587   if (x.type().IsCharacter()) {
588     CharacterMaxOrMin<false>(result, x, dim, source, line, mask, "MINVAL");
589   } else {
590     NumericMaxOrMin<false>(result, x, dim, source, line, mask, "MINVAL");
591   }
592 }
593 } // extern "C"
594 
595 // NORM2
596 
597 template <int KIND> class Norm2Accumulator {
598 public:
599   using Type = CppTypeFor<TypeCategory::Real, KIND>;
600   // Use at least double precision for accumulators
601   using AccumType = CppTypeFor<TypeCategory::Real, std::max(KIND, 8)>;
602   explicit Norm2Accumulator(const Descriptor &array) : array_{array} {}
603   void Reinitialize() { max_ = sum_ = 0; }
604   template <typename A> void GetResult(A *p, int /*zeroBasedDim*/ = -1) const {
605     // m * sqrt(1 + sum((others(:)/m)**2))
606     *p = static_cast<Type>(max_ * std::sqrt(1 + sum_));
607   }
608   bool Accumulate(Type x) {
609     auto absX{AccumType{std::abs(x)}};
610     if (!max_) {
611       max_ = x;
612     } else if (absX > max_) {
613       auto t{max_ / absX}; // < 1.0
614       auto tsq{t * t};
615       sum_ *= tsq; // scale sum to reflect change to the max
616       sum_ += tsq; // include a term for the previous max
617       max_ = absX;
618     } else { // absX <= max_
619       auto t{absX / max_};
620       sum_ += t * t;
621     }
622     return true;
623   }
624   template <typename A> bool AccumulateAt(const SubscriptValue at[]) {
625     return Accumulate(*array_.Element<A>(at));
626   }
627 
628 private:
629   const Descriptor &array_;
630   AccumType max_{0}; // value (m) with largest magnitude
631   AccumType sum_{0}; // sum((others(:)/m)**2)
632 };
633 
634 template <int KIND> struct Norm2Helper {
635   void operator()(Descriptor &result, const Descriptor &x, int dim,
636       const Descriptor *mask, Terminator &terminator) const {
637     DoMaxMinNorm2<TypeCategory::Real, KIND, Norm2Accumulator<KIND>>(
638         result, x, dim, mask, "NORM2", terminator);
639   }
640 };
641 
642 extern "C" {
643 // TODO: REAL(2 & 3)
644 CppTypeFor<TypeCategory::Real, 4> RTNAME(Norm2_4)(const Descriptor &x,
645     const char *source, int line, int dim, const Descriptor *mask) {
646   return GetTotalReduction<TypeCategory::Real, 4>(
647       x, source, line, dim, mask, Norm2Accumulator<4>{x}, "NORM2");
648 }
649 CppTypeFor<TypeCategory::Real, 8> RTNAME(Norm2_8)(const Descriptor &x,
650     const char *source, int line, int dim, const Descriptor *mask) {
651   return GetTotalReduction<TypeCategory::Real, 8>(
652       x, source, line, dim, mask, Norm2Accumulator<8>{x}, "NORM2");
653 }
654 #if LONG_DOUBLE == 80
655 CppTypeFor<TypeCategory::Real, 10> RTNAME(Norm2_10)(const Descriptor &x,
656     const char *source, int line, int dim, const Descriptor *mask) {
657   return GetTotalReduction<TypeCategory::Real, 10>(
658       x, source, line, dim, mask, Norm2Accumulator<10>{x}, "NORM2");
659 }
660 #elif LONG_DOUBLE == 128
661 CppTypeFor<TypeCategory::Real, 16> RTNAME(Norm2_16)(const Descriptor &x,
662     const char *source, int line, int dim, const Descriptor *mask) {
663   return GetTotalReduction<TypeCategory::Real, 16>(
664       x, source, line, dim, mask, Norm2Accumulator<16>{x}, "NORM2");
665 }
666 #endif
667 
668 void RTNAME(Norm2Dim)(Descriptor &result, const Descriptor &x, int dim,
669     const char *source, int line, const Descriptor *mask) {
670   Terminator terminator{source, line};
671   auto type{x.type().GetCategoryAndKind()};
672   RUNTIME_CHECK(terminator, type);
673   if (type->first == TypeCategory::Real) {
674     ApplyFloatingPointKind<Norm2Helper, void>(
675         type->second, terminator, result, x, dim, mask, terminator);
676   } else {
677     terminator.Crash("NORM2: bad type code %d", x.type().raw());
678   }
679 }
680 } // extern "C"
681 } // namespace Fortran::runtime
682