1 /* Implementation of the MAXLOC intrinsic
2 Copyright (C) 2002-2022 Free Software Foundation, Inc.
3 Contributed by Paul Brook <paul@nowt.org>
4
5 This file is part of the GNU Fortran runtime library (libgfortran).
6
7 Libgfortran is free software; you can redistribute it and/or
8 modify it under the terms of the GNU General Public
9 License as published by the Free Software Foundation; either
10 version 3 of the License, or (at your option) any later version.
11
12 Libgfortran is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 Under Section 7 of GPL version 3, you are granted additional
18 permissions described in the GCC Runtime Library Exception, version
19 3.1, as published by the Free Software Foundation.
20
21 You should have received a copy of the GNU General Public License and
22 a copy of the GCC Runtime Library Exception along with this program;
23 see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
24 <http://www.gnu.org/licenses/>. */
25
26 #include "libgfortran.h"
27 #include <assert.h>
28
29
30 #if defined (HAVE_GFC_REAL_8) && defined (HAVE_GFC_INTEGER_4)
31
32 #define HAVE_BACK_ARG 1
33
34
35 extern void maxloc1_4_r8 (gfc_array_i4 * const restrict,
36 gfc_array_r8 * const restrict, const index_type * const restrict, GFC_LOGICAL_4 back);
37 export_proto(maxloc1_4_r8);
38
39 void
maxloc1_4_r8(gfc_array_i4 * const restrict retarray,gfc_array_r8 * const restrict array,const index_type * const restrict pdim,GFC_LOGICAL_4 back)40 maxloc1_4_r8 (gfc_array_i4 * const restrict retarray,
41 gfc_array_r8 * const restrict array,
42 const index_type * const restrict pdim, GFC_LOGICAL_4 back)
43 {
44 index_type count[GFC_MAX_DIMENSIONS];
45 index_type extent[GFC_MAX_DIMENSIONS];
46 index_type sstride[GFC_MAX_DIMENSIONS];
47 index_type dstride[GFC_MAX_DIMENSIONS];
48 const GFC_REAL_8 * restrict base;
49 GFC_INTEGER_4 * restrict dest;
50 index_type rank;
51 index_type n;
52 index_type len;
53 index_type delta;
54 index_type dim;
55 int continue_loop;
56
57 /* Make dim zero based to avoid confusion. */
58 rank = GFC_DESCRIPTOR_RANK (array) - 1;
59 dim = (*pdim) - 1;
60
61 if (unlikely (dim < 0 || dim > rank))
62 {
63 runtime_error ("Dim argument incorrect in MAXLOC intrinsic: "
64 "is %ld, should be between 1 and %ld",
65 (long int) dim + 1, (long int) rank + 1);
66 }
67
68 len = GFC_DESCRIPTOR_EXTENT(array,dim);
69 if (len < 0)
70 len = 0;
71 delta = GFC_DESCRIPTOR_STRIDE(array,dim);
72
73 for (n = 0; n < dim; n++)
74 {
75 sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n);
76 extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
77
78 if (extent[n] < 0)
79 extent[n] = 0;
80 }
81 for (n = dim; n < rank; n++)
82 {
83 sstride[n] = GFC_DESCRIPTOR_STRIDE(array, n + 1);
84 extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1);
85
86 if (extent[n] < 0)
87 extent[n] = 0;
88 }
89
90 if (retarray->base_addr == NULL)
91 {
92 size_t alloc_size, str;
93
94 for (n = 0; n < rank; n++)
95 {
96 if (n == 0)
97 str = 1;
98 else
99 str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
100
101 GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);
102
103 }
104
105 retarray->offset = 0;
106 retarray->dtype.rank = rank;
107
108 alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1];
109
110 retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_4));
111 if (alloc_size == 0)
112 {
113 /* Make sure we have a zero-sized array. */
114 GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
115 return;
116
117 }
118 }
119 else
120 {
121 if (rank != GFC_DESCRIPTOR_RANK (retarray))
122 runtime_error ("rank of return array incorrect in"
123 " MAXLOC intrinsic: is %ld, should be %ld",
124 (long int) (GFC_DESCRIPTOR_RANK (retarray)),
125 (long int) rank);
126
127 if (unlikely (compile_options.bounds_check))
128 bounds_ifunction_return ((array_t *) retarray, extent,
129 "return value", "MAXLOC");
130 }
131
132 for (n = 0; n < rank; n++)
133 {
134 count[n] = 0;
135 dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
136 if (extent[n] <= 0)
137 return;
138 }
139
140 base = array->base_addr;
141 dest = retarray->base_addr;
142
143 continue_loop = 1;
144 while (continue_loop)
145 {
146 const GFC_REAL_8 * restrict src;
147 GFC_INTEGER_4 result;
148 src = base;
149 {
150
151 GFC_REAL_8 maxval;
152 #if defined (GFC_REAL_8_INFINITY)
153 maxval = -GFC_REAL_8_INFINITY;
154 #else
155 maxval = -GFC_REAL_8_HUGE;
156 #endif
157 result = 1;
158 if (len <= 0)
159 *dest = 0;
160 else
161 {
162 #if ! defined HAVE_BACK_ARG
163 for (n = 0; n < len; n++, src += delta)
164 {
165 #endif
166
167 #if defined (GFC_REAL_8_QUIET_NAN)
168 for (n = 0; n < len; n++, src += delta)
169 {
170 if (*src >= maxval)
171 {
172 maxval = *src;
173 result = (GFC_INTEGER_4)n + 1;
174 break;
175 }
176 }
177 #else
178 n = 0;
179 #endif
180 for (; n < len; n++, src += delta)
181 {
182 if (back ? *src >= maxval : *src > maxval)
183 {
184 maxval = *src;
185 result = (GFC_INTEGER_4)n + 1;
186 }
187 }
188
189 *dest = result;
190 }
191 }
192 /* Advance to the next element. */
193 count[0]++;
194 base += sstride[0];
195 dest += dstride[0];
196 n = 0;
197 while (count[n] == extent[n])
198 {
199 /* When we get to the end of a dimension, reset it and increment
200 the next dimension. */
201 count[n] = 0;
202 /* We could precalculate these products, but this is a less
203 frequently used path so probably not worth it. */
204 base -= sstride[n] * extent[n];
205 dest -= dstride[n] * extent[n];
206 n++;
207 if (n >= rank)
208 {
209 /* Break out of the loop. */
210 continue_loop = 0;
211 break;
212 }
213 else
214 {
215 count[n]++;
216 base += sstride[n];
217 dest += dstride[n];
218 }
219 }
220 }
221 }
222
223
224 extern void mmaxloc1_4_r8 (gfc_array_i4 * const restrict,
225 gfc_array_r8 * const restrict, const index_type * const restrict,
226 gfc_array_l1 * const restrict, GFC_LOGICAL_4 back);
227 export_proto(mmaxloc1_4_r8);
228
229 void
230 mmaxloc1_4_r8 (gfc_array_i4 * const restrict retarray,
231 gfc_array_r8 * const restrict array,
232 const index_type * const restrict pdim,
233 gfc_array_l1 * const restrict mask, GFC_LOGICAL_4 back)
234 {
235 index_type count[GFC_MAX_DIMENSIONS];
236 index_type extent[GFC_MAX_DIMENSIONS];
237 index_type sstride[GFC_MAX_DIMENSIONS];
238 index_type dstride[GFC_MAX_DIMENSIONS];
239 index_type mstride[GFC_MAX_DIMENSIONS];
240 GFC_INTEGER_4 * restrict dest;
241 const GFC_REAL_8 * restrict base;
242 const GFC_LOGICAL_1 * restrict mbase;
243 index_type rank;
244 index_type dim;
245 index_type n;
246 index_type len;
247 index_type delta;
248 index_type mdelta;
249 int mask_kind;
250
251 if (mask == NULL)
252 {
253 #ifdef HAVE_BACK_ARG
254 maxloc1_4_r8 (retarray, array, pdim, back);
255 #else
256 maxloc1_4_r8 (retarray, array, pdim);
257 #endif
258 return;
259 }
260
261 dim = (*pdim) - 1;
262 rank = GFC_DESCRIPTOR_RANK (array) - 1;
263
264
265 if (unlikely (dim < 0 || dim > rank))
266 {
267 runtime_error ("Dim argument incorrect in MAXLOC intrinsic: "
268 "is %ld, should be between 1 and %ld",
269 (long int) dim + 1, (long int) rank + 1);
270 }
271
272 len = GFC_DESCRIPTOR_EXTENT(array,dim);
273 if (len <= 0)
274 return;
275
276 mbase = mask->base_addr;
277
278 mask_kind = GFC_DESCRIPTOR_SIZE (mask);
279
280 if (mask_kind == 1 || mask_kind == 2 || mask_kind == 4 || mask_kind == 8
281 #ifdef HAVE_GFC_LOGICAL_16
282 || mask_kind == 16
283 #endif
284 )
285 mbase = GFOR_POINTER_TO_L1 (mbase, mask_kind);
286 else
287 runtime_error ("Funny sized logical array");
288
289 delta = GFC_DESCRIPTOR_STRIDE(array,dim);
290 mdelta = GFC_DESCRIPTOR_STRIDE_BYTES(mask,dim);
291
292 for (n = 0; n < dim; n++)
293 {
294 sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n);
295 mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask,n);
296 extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
297
298 if (extent[n] < 0)
299 extent[n] = 0;
300
301 }
302 for (n = dim; n < rank; n++)
303 {
304 sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n + 1);
305 mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask, n + 1);
306 extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1);
307
308 if (extent[n] < 0)
309 extent[n] = 0;
310 }
311
312 if (retarray->base_addr == NULL)
313 {
314 size_t alloc_size, str;
315
316 for (n = 0; n < rank; n++)
317 {
318 if (n == 0)
319 str = 1;
320 else
321 str= GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
322
323 GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);
324
325 }
326
327 alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1];
328
329 retarray->offset = 0;
330 retarray->dtype.rank = rank;
331
332 if (alloc_size == 0)
333 {
334 /* Make sure we have a zero-sized array. */
335 GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
336 return;
337 }
338 else
339 retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_4));
340
341 }
342 else
343 {
344 if (rank != GFC_DESCRIPTOR_RANK (retarray))
345 runtime_error ("rank of return array incorrect in MAXLOC intrinsic");
346
347 if (unlikely (compile_options.bounds_check))
348 {
349 bounds_ifunction_return ((array_t *) retarray, extent,
350 "return value", "MAXLOC");
351 bounds_equal_extents ((array_t *) mask, (array_t *) array,
352 "MASK argument", "MAXLOC");
353 }
354 }
355
356 for (n = 0; n < rank; n++)
357 {
358 count[n] = 0;
359 dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
360 if (extent[n] <= 0)
361 return;
362 }
363
364 dest = retarray->base_addr;
365 base = array->base_addr;
366
367 while (base)
368 {
369 const GFC_REAL_8 * restrict src;
370 const GFC_LOGICAL_1 * restrict msrc;
371 GFC_INTEGER_4 result;
372 src = base;
373 msrc = mbase;
374 {
375
376 GFC_REAL_8 maxval;
377 #if defined (GFC_REAL_8_INFINITY)
378 maxval = -GFC_REAL_8_INFINITY;
379 #else
380 maxval = -GFC_REAL_8_HUGE;
381 #endif
382 #if defined (GFC_REAL_8_QUIET_NAN)
383 GFC_INTEGER_4 result2 = 0;
384 #endif
385 result = 0;
386 for (n = 0; n < len; n++, src += delta, msrc += mdelta)
387 {
388
389 if (*msrc)
390 {
391 #if defined (GFC_REAL_8_QUIET_NAN)
392 if (!result2)
393 result2 = (GFC_INTEGER_4)n + 1;
394 if (*src >= maxval)
395 #endif
396 {
397 maxval = *src;
398 result = (GFC_INTEGER_4)n + 1;
399 break;
400 }
401 }
402 }
403 #if defined (GFC_REAL_8_QUIET_NAN)
404 if (unlikely (n >= len))
405 result = result2;
406 else
407 #endif
408 if (back)
409 for (; n < len; n++, src += delta, msrc += mdelta)
410 {
411 if (*msrc && unlikely (*src >= maxval))
412 {
413 maxval = *src;
414 result = (GFC_INTEGER_4)n + 1;
415 }
416 }
417 else
418 for (; n < len; n++, src += delta, msrc += mdelta)
419 {
420 if (*msrc && unlikely (*src > maxval))
421 {
422 maxval = *src;
423 result = (GFC_INTEGER_4)n + 1;
424 }
425 }
426 *dest = result;
427 }
428 /* Advance to the next element. */
429 count[0]++;
430 base += sstride[0];
431 mbase += mstride[0];
432 dest += dstride[0];
433 n = 0;
434 while (count[n] == extent[n])
435 {
436 /* When we get to the end of a dimension, reset it and increment
437 the next dimension. */
438 count[n] = 0;
439 /* We could precalculate these products, but this is a less
440 frequently used path so probably not worth it. */
441 base -= sstride[n] * extent[n];
442 mbase -= mstride[n] * extent[n];
443 dest -= dstride[n] * extent[n];
444 n++;
445 if (n >= rank)
446 {
447 /* Break out of the loop. */
448 base = NULL;
449 break;
450 }
451 else
452 {
453 count[n]++;
454 base += sstride[n];
455 mbase += mstride[n];
456 dest += dstride[n];
457 }
458 }
459 }
460 }
461
462
463 extern void smaxloc1_4_r8 (gfc_array_i4 * const restrict,
464 gfc_array_r8 * const restrict, const index_type * const restrict,
465 GFC_LOGICAL_4 *, GFC_LOGICAL_4 back);
466 export_proto(smaxloc1_4_r8);
467
468 void
469 smaxloc1_4_r8 (gfc_array_i4 * const restrict retarray,
470 gfc_array_r8 * const restrict array,
471 const index_type * const restrict pdim,
472 GFC_LOGICAL_4 * mask, GFC_LOGICAL_4 back)
473 {
474 index_type count[GFC_MAX_DIMENSIONS];
475 index_type extent[GFC_MAX_DIMENSIONS];
476 index_type dstride[GFC_MAX_DIMENSIONS];
477 GFC_INTEGER_4 * restrict dest;
478 index_type rank;
479 index_type n;
480 index_type dim;
481
482
483 if (mask == NULL || *mask)
484 {
485 #ifdef HAVE_BACK_ARG
486 maxloc1_4_r8 (retarray, array, pdim, back);
487 #else
488 maxloc1_4_r8 (retarray, array, pdim);
489 #endif
490 return;
491 }
492 /* Make dim zero based to avoid confusion. */
493 dim = (*pdim) - 1;
494 rank = GFC_DESCRIPTOR_RANK (array) - 1;
495
496 if (unlikely (dim < 0 || dim > rank))
497 {
498 runtime_error ("Dim argument incorrect in MAXLOC intrinsic: "
499 "is %ld, should be between 1 and %ld",
500 (long int) dim + 1, (long int) rank + 1);
501 }
502
503 for (n = 0; n < dim; n++)
504 {
505 extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
506
507 if (extent[n] <= 0)
508 extent[n] = 0;
509 }
510
511 for (n = dim; n < rank; n++)
512 {
513 extent[n] =
514 GFC_DESCRIPTOR_EXTENT(array,n + 1);
515
516 if (extent[n] <= 0)
517 extent[n] = 0;
518 }
519
520 if (retarray->base_addr == NULL)
521 {
522 size_t alloc_size, str;
523
524 for (n = 0; n < rank; n++)
525 {
526 if (n == 0)
527 str = 1;
528 else
529 str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
530
531 GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);
532
533 }
534
535 retarray->offset = 0;
536 retarray->dtype.rank = rank;
537
538 alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1];
539
540 if (alloc_size == 0)
541 {
542 /* Make sure we have a zero-sized array. */
543 GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
544 return;
545 }
546 else
547 retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_4));
548 }
549 else
550 {
551 if (rank != GFC_DESCRIPTOR_RANK (retarray))
552 runtime_error ("rank of return array incorrect in"
553 " MAXLOC intrinsic: is %ld, should be %ld",
554 (long int) (GFC_DESCRIPTOR_RANK (retarray)),
555 (long int) rank);
556
557 if (unlikely (compile_options.bounds_check))
558 {
559 for (n=0; n < rank; n++)
560 {
561 index_type ret_extent;
562
563 ret_extent = GFC_DESCRIPTOR_EXTENT(retarray,n);
564 if (extent[n] != ret_extent)
565 runtime_error ("Incorrect extent in return value of"
566 " MAXLOC intrinsic in dimension %ld:"
567 " is %ld, should be %ld", (long int) n + 1,
568 (long int) ret_extent, (long int) extent[n]);
569 }
570 }
571 }
572
573 for (n = 0; n < rank; n++)
574 {
575 count[n] = 0;
576 dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
577 }
578
579 dest = retarray->base_addr;
580
581 while(1)
582 {
583 *dest = 0;
584 count[0]++;
585 dest += dstride[0];
586 n = 0;
587 while (count[n] == extent[n])
588 {
589 /* When we get to the end of a dimension, reset it and increment
590 the next dimension. */
591 count[n] = 0;
592 /* We could precalculate these products, but this is a less
593 frequently used path so probably not worth it. */
594 dest -= dstride[n] * extent[n];
595 n++;
596 if (n >= rank)
597 return;
598 else
599 {
600 count[n]++;
601 dest += dstride[n];
602 }
603 }
604 }
605 }
606
607 #endif
608