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