197403Sobrien/* Implementation of the PRODUCT intrinsic
297403Sobrien   Copyright (C) 2002-2020 Free Software Foundation, Inc.
397403Sobrien   Contributed by Paul Brook <paul@nowt.org>
497403Sobrien
597403SobrienThis file is part of the GNU Fortran runtime library (libgfortran).
697403Sobrien
797403SobrienLibgfortran is free software; you can redistribute it and/or
897403Sobrienmodify it under the terms of the GNU General Public
997403SobrienLicense as published by the Free Software Foundation; either
1097403Sobrienversion 3 of the License, or (at your option) any later version.
1197403Sobrien
1297403SobrienLibgfortran is distributed in the hope that it will be useful,
1397403Sobrienbut WITHOUT ANY WARRANTY; without even the implied warranty of
1497403SobrienMERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
1597403SobrienGNU General Public License for more details.
1697403Sobrien
1797403SobrienUnder Section 7 of GPL version 3, you are granted additional
1897403Sobrienpermissions described in the GCC Runtime Library Exception, version
1997403Sobrien3.1, as published by the Free Software Foundation.
2097403Sobrien
2197403SobrienYou should have received a copy of the GNU General Public License and
2297403Sobriena copy of the GCC Runtime Library Exception along with this program;
2397403Sobriensee the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
2497403Sobrien<http://www.gnu.org/licenses/>.  */
2597403Sobrien
2697403Sobrien#include "libgfortran.h"
2797403Sobrien
2897403Sobrien
2997403Sobrien#if defined (HAVE_GFC_INTEGER_4) && defined (HAVE_GFC_INTEGER_4)
3097403Sobrien
3197403Sobrien
3297403Sobrienextern void product_i4 (gfc_array_i4 * const restrict,
3397403Sobrien	gfc_array_i4 * const restrict, const index_type * const restrict);
3497403Sobrienexport_proto(product_i4);
3597403Sobrien
3697403Sobrienvoid
3797403Sobrienproduct_i4 (gfc_array_i4 * const restrict retarray,
3897403Sobrien	gfc_array_i4 * const restrict array,
3997403Sobrien	const index_type * const restrict pdim)
4097403Sobrien{
4197403Sobrien  index_type count[GFC_MAX_DIMENSIONS];
4297403Sobrien  index_type extent[GFC_MAX_DIMENSIONS];
4397403Sobrien  index_type sstride[GFC_MAX_DIMENSIONS];
4497403Sobrien  index_type dstride[GFC_MAX_DIMENSIONS];
4597403Sobrien  const GFC_INTEGER_4 * restrict base;
4697403Sobrien  GFC_INTEGER_4 * restrict dest;
4797403Sobrien  index_type rank;
4897403Sobrien  index_type n;
4997403Sobrien  index_type len;
5097403Sobrien  index_type delta;
5197403Sobrien  index_type dim;
5297403Sobrien  int continue_loop;
5397403Sobrien
5497403Sobrien  /* Make dim zero based to avoid confusion.  */
5597403Sobrien  rank = GFC_DESCRIPTOR_RANK (array) - 1;
5697403Sobrien  dim = (*pdim) - 1;
5797403Sobrien
5897403Sobrien  if (unlikely (dim < 0 || dim > rank))
5997403Sobrien    {
6097403Sobrien      runtime_error ("Dim argument incorrect in PRODUCT intrinsic: "
6197403Sobrien 		     "is %ld, should be between 1 and %ld",
6297403Sobrien		     (long int) dim + 1, (long int) rank + 1);
6397403Sobrien    }
6497403Sobrien
6597403Sobrien  len = GFC_DESCRIPTOR_EXTENT(array,dim);
6697403Sobrien  if (len < 0)
6797403Sobrien    len = 0;
6897403Sobrien  delta = GFC_DESCRIPTOR_STRIDE(array,dim);
6997403Sobrien
7097403Sobrien  for (n = 0; n < dim; n++)
7197403Sobrien    {
7297403Sobrien      sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n);
7397403Sobrien      extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
7497403Sobrien
7597403Sobrien      if (extent[n] < 0)
7697403Sobrien	extent[n] = 0;
7797403Sobrien    }
7897403Sobrien  for (n = dim; n < rank; n++)
79110614Skan    {
8097403Sobrien      sstride[n] = GFC_DESCRIPTOR_STRIDE(array, n + 1);
8197403Sobrien      extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1);
8297403Sobrien
8397403Sobrien      if (extent[n] < 0)
8497403Sobrien	extent[n] = 0;
8597403Sobrien    }
8697403Sobrien
8797403Sobrien  if (retarray->base_addr == NULL)
8897403Sobrien    {
8997403Sobrien      size_t alloc_size, str;
9097403Sobrien
9197403Sobrien      for (n = 0; n < rank; n++)
9297403Sobrien	{
9397403Sobrien	  if (n == 0)
9497403Sobrien	    str = 1;
9597403Sobrien	  else
9697403Sobrien	    str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
9797403Sobrien
9897403Sobrien	  GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);
9997403Sobrien
10097403Sobrien	}
10197403Sobrien
10297403Sobrien      retarray->offset = 0;
10397403Sobrien      retarray->dtype.rank = rank;
10497403Sobrien
10597403Sobrien      alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1];
10697403Sobrien
10797403Sobrien      retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_4));
10897403Sobrien      if (alloc_size == 0)
10997403Sobrien	{
11097403Sobrien	  /* Make sure we have a zero-sized array.  */
11197403Sobrien	  GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
11297403Sobrien	  return;
11397403Sobrien
114102782Skan	}
11597403Sobrien    }
11697403Sobrien  else
11797403Sobrien    {
11897403Sobrien      if (rank != GFC_DESCRIPTOR_RANK (retarray))
11997403Sobrien	runtime_error ("rank of return array incorrect in"
12097403Sobrien		       " PRODUCT intrinsic: is %ld, should be %ld",
12197403Sobrien		       (long int) (GFC_DESCRIPTOR_RANK (retarray)),
12297403Sobrien		       (long int) rank);
12397403Sobrien
12497403Sobrien      if (unlikely (compile_options.bounds_check))
12597403Sobrien	bounds_ifunction_return ((array_t *) retarray, extent,
12697403Sobrien				 "return value", "PRODUCT");
12797403Sobrien    }
12897403Sobrien
12997403Sobrien  for (n = 0; n < rank; n++)
13097403Sobrien    {
13197403Sobrien      count[n] = 0;
13297403Sobrien      dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
13397403Sobrien      if (extent[n] <= 0)
13497403Sobrien	return;
13597403Sobrien    }
13697403Sobrien
13797403Sobrien  base = array->base_addr;
13897403Sobrien  dest = retarray->base_addr;
13997403Sobrien
14097403Sobrien  continue_loop = 1;
14197403Sobrien  while (continue_loop)
14297403Sobrien    {
14397403Sobrien      const GFC_INTEGER_4 * restrict src;
14497403Sobrien      GFC_INTEGER_4 result;
14597403Sobrien      src = base;
14697403Sobrien      {
14797403Sobrien
14897403Sobrien  result = 1;
14997403Sobrien	if (len <= 0)
15097403Sobrien	  *dest = 1;
15197403Sobrien	else
15297403Sobrien	  {
153103447Skan#if ! defined HAVE_BACK_ARG
154103447Skan	    for (n = 0; n < len; n++, src += delta)
155103447Skan	      {
15697403Sobrien#endif
15797403Sobrien
15897403Sobrien  result *= *src;
15997403Sobrien	      }
16097403Sobrien
16197403Sobrien	    *dest = result;
16297403Sobrien	  }
16397403Sobrien      }
16497403Sobrien      /* Advance to the next element.  */
16597403Sobrien      count[0]++;
16697403Sobrien      base += sstride[0];
16797403Sobrien      dest += dstride[0];
16897403Sobrien      n = 0;
16997403Sobrien      while (count[n] == extent[n])
17097403Sobrien	{
171103447Skan	  /* When we get to the end of a dimension, reset it and increment
17297403Sobrien	     the next dimension.  */
17397403Sobrien	  count[n] = 0;
17497403Sobrien	  /* We could precalculate these products, but this is a less
17597403Sobrien	     frequently used path so probably not worth it.  */
17697403Sobrien	  base -= sstride[n] * extent[n];
17797403Sobrien	  dest -= dstride[n] * extent[n];
17897403Sobrien	  n++;
17997403Sobrien	  if (n >= rank)
18097403Sobrien	    {
18197403Sobrien	      /* Break out of the loop.  */
18297403Sobrien	      continue_loop = 0;
18397403Sobrien	      break;
18497403Sobrien	    }
18597403Sobrien	  else
18697403Sobrien	    {
18797403Sobrien	      count[n]++;
18897403Sobrien	      base += sstride[n];
18997403Sobrien	      dest += dstride[n];
19097403Sobrien	    }
19197403Sobrien	}
19297403Sobrien    }
19397403Sobrien}
19497403Sobrien
19597403Sobrien
19697403Sobrienextern void mproduct_i4 (gfc_array_i4 * const restrict,
19797403Sobrien	gfc_array_i4 * const restrict, const index_type * const restrict,
19897403Sobrien	gfc_array_l1 * const restrict);
19997403Sobrienexport_proto(mproduct_i4);
20097403Sobrien
20197403Sobrienvoid
20297403Sobrienmproduct_i4 (gfc_array_i4 * const restrict retarray,
20397403Sobrien	gfc_array_i4 * const restrict array,
20497403Sobrien	const index_type * const restrict pdim,
20597403Sobrien	gfc_array_l1 * const restrict mask)
20697403Sobrien{
20797403Sobrien  index_type count[GFC_MAX_DIMENSIONS];
20897403Sobrien  index_type extent[GFC_MAX_DIMENSIONS];
20997403Sobrien  index_type sstride[GFC_MAX_DIMENSIONS];
21097403Sobrien  index_type dstride[GFC_MAX_DIMENSIONS];
21197403Sobrien  index_type mstride[GFC_MAX_DIMENSIONS];
21297403Sobrien  GFC_INTEGER_4 * restrict dest;
21397403Sobrien  const GFC_INTEGER_4 * restrict base;
21497403Sobrien  const GFC_LOGICAL_1 * restrict mbase;
21597403Sobrien  index_type rank;
21697403Sobrien  index_type dim;
21797403Sobrien  index_type n;
21897403Sobrien  index_type len;
21997403Sobrien  index_type delta;
22097403Sobrien  index_type mdelta;
22197403Sobrien  int mask_kind;
22297403Sobrien
22397403Sobrien  if (mask == NULL)
22497403Sobrien    {
22597403Sobrien#ifdef HAVE_BACK_ARG
22697403Sobrien      product_i4 (retarray, array, pdim, back);
22797403Sobrien#else
22897403Sobrien      product_i4 (retarray, array, pdim);
22997403Sobrien#endif
23097403Sobrien      return;
23197403Sobrien    }
23297403Sobrien
23397403Sobrien  dim = (*pdim) - 1;
23497403Sobrien  rank = GFC_DESCRIPTOR_RANK (array) - 1;
23597403Sobrien
23697403Sobrien
23797403Sobrien  if (unlikely (dim < 0 || dim > rank))
23897403Sobrien    {
23997403Sobrien      runtime_error ("Dim argument incorrect in PRODUCT intrinsic: "
24097403Sobrien 		     "is %ld, should be between 1 and %ld",
24197403Sobrien		     (long int) dim + 1, (long int) rank + 1);
24297403Sobrien    }
24397403Sobrien
24497403Sobrien  len = GFC_DESCRIPTOR_EXTENT(array,dim);
24597403Sobrien  if (len <= 0)
24697403Sobrien    return;
24797403Sobrien
24897403Sobrien  mbase = mask->base_addr;
24997403Sobrien
25097403Sobrien  mask_kind = GFC_DESCRIPTOR_SIZE (mask);
25197403Sobrien
25297403Sobrien  if (mask_kind == 1 || mask_kind == 2 || mask_kind == 4 || mask_kind == 8
25397403Sobrien#ifdef HAVE_GFC_LOGICAL_16
25497403Sobrien      || mask_kind == 16
25597403Sobrien#endif
25697403Sobrien      )
25797403Sobrien    mbase = GFOR_POINTER_TO_L1 (mbase, mask_kind);
25897403Sobrien  else
25997403Sobrien    runtime_error ("Funny sized logical array");
26097403Sobrien
26197403Sobrien  delta = GFC_DESCRIPTOR_STRIDE(array,dim);
26297403Sobrien  mdelta = GFC_DESCRIPTOR_STRIDE_BYTES(mask,dim);
26397403Sobrien
26497403Sobrien  for (n = 0; n < dim; n++)
26597403Sobrien    {
26697403Sobrien      sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n);
26797403Sobrien      mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask,n);
26897403Sobrien      extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
26997403Sobrien
27097403Sobrien      if (extent[n] < 0)
27197403Sobrien	extent[n] = 0;
27297403Sobrien
27397403Sobrien    }
27497403Sobrien  for (n = dim; n < rank; n++)
27597403Sobrien    {
27697403Sobrien      sstride[n] = GFC_DESCRIPTOR_STRIDE(array,n + 1);
27797403Sobrien      mstride[n] = GFC_DESCRIPTOR_STRIDE_BYTES(mask, n + 1);
27897403Sobrien      extent[n] = GFC_DESCRIPTOR_EXTENT(array, n + 1);
27997403Sobrien
28097403Sobrien      if (extent[n] < 0)
28197403Sobrien	extent[n] = 0;
28297403Sobrien    }
28397403Sobrien
28497403Sobrien  if (retarray->base_addr == NULL)
28597403Sobrien    {
28697403Sobrien      size_t alloc_size, str;
28797403Sobrien
28897403Sobrien      for (n = 0; n < rank; n++)
28997403Sobrien	{
29097403Sobrien	  if (n == 0)
29197403Sobrien	    str = 1;
29297403Sobrien	  else
29397403Sobrien	    str= GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
29497403Sobrien
29597403Sobrien	  GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);
29697403Sobrien
29797403Sobrien	}
29897403Sobrien
29997403Sobrien      alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1];
30097403Sobrien
30197403Sobrien      retarray->offset = 0;
30297403Sobrien      retarray->dtype.rank = rank;
30397403Sobrien
30497403Sobrien      if (alloc_size == 0)
30597403Sobrien	{
30697403Sobrien	  /* Make sure we have a zero-sized array.  */
30797403Sobrien	  GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
30897403Sobrien	  return;
30997403Sobrien	}
31097403Sobrien      else
31197403Sobrien	retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_4));
31297403Sobrien
31397403Sobrien    }
31497403Sobrien  else
31597403Sobrien    {
31697403Sobrien      if (rank != GFC_DESCRIPTOR_RANK (retarray))
31797403Sobrien	runtime_error ("rank of return array incorrect in PRODUCT intrinsic");
31897403Sobrien
31997403Sobrien      if (unlikely (compile_options.bounds_check))
32097403Sobrien	{
32197403Sobrien	  bounds_ifunction_return ((array_t *) retarray, extent,
32297403Sobrien				   "return value", "PRODUCT");
32397403Sobrien	  bounds_equal_extents ((array_t *) mask, (array_t *) array,
32497403Sobrien	  			"MASK argument", "PRODUCT");
32597403Sobrien	}
32697403Sobrien    }
32797403Sobrien
32897403Sobrien  for (n = 0; n < rank; n++)
32997403Sobrien    {
33097403Sobrien      count[n] = 0;
33197403Sobrien      dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
33297403Sobrien      if (extent[n] <= 0)
33397403Sobrien	return;
33497403Sobrien    }
33597403Sobrien
33697403Sobrien  dest = retarray->base_addr;
33797403Sobrien  base = array->base_addr;
33897403Sobrien
33997403Sobrien  while (base)
34097403Sobrien    {
34197403Sobrien      const GFC_INTEGER_4 * restrict src;
34297403Sobrien      const GFC_LOGICAL_1 * restrict msrc;
34397403Sobrien      GFC_INTEGER_4 result;
34497403Sobrien      src = base;
34597403Sobrien      msrc = mbase;
34697403Sobrien      {
34797403Sobrien
34897403Sobrien  result = 1;
34997403Sobrien	for (n = 0; n < len; n++, src += delta, msrc += mdelta)
35097403Sobrien	  {
35197403Sobrien
35297403Sobrien  if (*msrc)
35397403Sobrien    result *= *src;
35497403Sobrien	  }
35597403Sobrien	*dest = result;
35697403Sobrien      }
35797403Sobrien      /* Advance to the next element.  */
35897403Sobrien      count[0]++;
35997403Sobrien      base += sstride[0];
36097403Sobrien      mbase += mstride[0];
36197403Sobrien      dest += dstride[0];
36297403Sobrien      n = 0;
36397403Sobrien      while (count[n] == extent[n])
36497403Sobrien	{
36597403Sobrien	  /* When we get to the end of a dimension, reset it and increment
36697403Sobrien	     the next dimension.  */
36797403Sobrien	  count[n] = 0;
36897403Sobrien	  /* We could precalculate these products, but this is a less
36997403Sobrien	     frequently used path so probably not worth it.  */
37097403Sobrien	  base -= sstride[n] * extent[n];
37197403Sobrien	  mbase -= mstride[n] * extent[n];
37297403Sobrien	  dest -= dstride[n] * extent[n];
373	  n++;
374	  if (n >= rank)
375	    {
376	      /* Break out of the loop.  */
377	      base = NULL;
378	      break;
379	    }
380	  else
381	    {
382	      count[n]++;
383	      base += sstride[n];
384	      mbase += mstride[n];
385	      dest += dstride[n];
386	    }
387	}
388    }
389}
390
391
392extern void sproduct_i4 (gfc_array_i4 * const restrict,
393	gfc_array_i4 * const restrict, const index_type * const restrict,
394	GFC_LOGICAL_4 *);
395export_proto(sproduct_i4);
396
397void
398sproduct_i4 (gfc_array_i4 * const restrict retarray,
399	gfc_array_i4 * const restrict array,
400	const index_type * const restrict pdim,
401	GFC_LOGICAL_4 * mask)
402{
403  index_type count[GFC_MAX_DIMENSIONS];
404  index_type extent[GFC_MAX_DIMENSIONS];
405  index_type dstride[GFC_MAX_DIMENSIONS];
406  GFC_INTEGER_4 * restrict dest;
407  index_type rank;
408  index_type n;
409  index_type dim;
410
411
412  if (mask == NULL || *mask)
413    {
414#ifdef HAVE_BACK_ARG
415      product_i4 (retarray, array, pdim, back);
416#else
417      product_i4 (retarray, array, pdim);
418#endif
419      return;
420    }
421  /* Make dim zero based to avoid confusion.  */
422  dim = (*pdim) - 1;
423  rank = GFC_DESCRIPTOR_RANK (array) - 1;
424
425  if (unlikely (dim < 0 || dim > rank))
426    {
427      runtime_error ("Dim argument incorrect in PRODUCT intrinsic: "
428 		     "is %ld, should be between 1 and %ld",
429		     (long int) dim + 1, (long int) rank + 1);
430    }
431
432  for (n = 0; n < dim; n++)
433    {
434      extent[n] = GFC_DESCRIPTOR_EXTENT(array,n);
435
436      if (extent[n] <= 0)
437	extent[n] = 0;
438    }
439
440  for (n = dim; n < rank; n++)
441    {
442      extent[n] =
443	GFC_DESCRIPTOR_EXTENT(array,n + 1);
444
445      if (extent[n] <= 0)
446	extent[n] = 0;
447    }
448
449  if (retarray->base_addr == NULL)
450    {
451      size_t alloc_size, str;
452
453      for (n = 0; n < rank; n++)
454	{
455	  if (n == 0)
456	    str = 1;
457	  else
458	    str = GFC_DESCRIPTOR_STRIDE(retarray,n-1) * extent[n-1];
459
460	  GFC_DIMENSION_SET(retarray->dim[n], 0, extent[n] - 1, str);
461
462	}
463
464      retarray->offset = 0;
465      retarray->dtype.rank = rank;
466
467      alloc_size = GFC_DESCRIPTOR_STRIDE(retarray,rank-1) * extent[rank-1];
468
469      if (alloc_size == 0)
470	{
471	  /* Make sure we have a zero-sized array.  */
472	  GFC_DIMENSION_SET(retarray->dim[0], 0, -1, 1);
473	  return;
474	}
475      else
476	retarray->base_addr = xmallocarray (alloc_size, sizeof (GFC_INTEGER_4));
477    }
478  else
479    {
480      if (rank != GFC_DESCRIPTOR_RANK (retarray))
481	runtime_error ("rank of return array incorrect in"
482		       " PRODUCT intrinsic: is %ld, should be %ld",
483		       (long int) (GFC_DESCRIPTOR_RANK (retarray)),
484		       (long int) rank);
485
486      if (unlikely (compile_options.bounds_check))
487	{
488	  for (n=0; n < rank; n++)
489	    {
490	      index_type ret_extent;
491
492	      ret_extent = GFC_DESCRIPTOR_EXTENT(retarray,n);
493	      if (extent[n] != ret_extent)
494		runtime_error ("Incorrect extent in return value of"
495			       " PRODUCT intrinsic in dimension %ld:"
496			       " is %ld, should be %ld", (long int) n + 1,
497			       (long int) ret_extent, (long int) extent[n]);
498	    }
499	}
500    }
501
502  for (n = 0; n < rank; n++)
503    {
504      count[n] = 0;
505      dstride[n] = GFC_DESCRIPTOR_STRIDE(retarray,n);
506    }
507
508  dest = retarray->base_addr;
509
510  while(1)
511    {
512      *dest = 1;
513      count[0]++;
514      dest += dstride[0];
515      n = 0;
516      while (count[n] == extent[n])
517	{
518	  /* When we get to the end of a dimension, reset it and increment
519	     the next dimension.  */
520	  count[n] = 0;
521	  /* We could precalculate these products, but this is a less
522	     frequently used path so probably not worth it.  */
523	  dest -= dstride[n] * extent[n];
524	  n++;
525	  if (n >= rank)
526	    return;
527	  else
528	    {
529	      count[n]++;
530	      dest += dstride[n];
531	    }
532      	}
533    }
534}
535
536#endif
537