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