[thirdparty/gcc.git] / libgfortran / generated / spread_r8.c

/* Special implementation of the SPREAD intrinsic
   Copyright 2008, 2009 Free Software Foundation, Inc.
   Contributed by Thomas Koenig <tkoenig@gcc.gnu.org>, based on
   spread_generic.c written by Paul Brook <paul@nowt.org>

This file is part of the GNU Fortran 95 runtime library (libgfortran).

Libgfortran is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public
License as published by the Free Software Foundation; either
version 3 of the License, or (at your option) any later version.

Ligbfortran is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.

Under Section 7 of GPL version 3, you are granted additional
permissions described in the GCC Runtime Library Exception, version
3.1, as published by the Free Software Foundation.

You should have received a copy of the GNU General Public License and
a copy of the GCC Runtime Library Exception along with this program;
see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
<http://www.gnu.org/licenses/>.  */

#include "libgfortran.h"
#include <stdlib.h>
#include <assert.h>
#include <string.h>


#if defined (HAVE_GFC_REAL_8)

void
spread_r8 (gfc_array_r8 *ret, const gfc_array_r8 *source,
		 const index_type along, const index_type pncopies)
{
  /* r.* indicates the return array.  */
  index_type rstride[GFC_MAX_DIMENSIONS];
  index_type rstride0;
  index_type rdelta = 0;
  index_type rrank;
  index_type rs;
  GFC_REAL_8 *rptr;
  GFC_REAL_8 * restrict dest;
  /* s.* indicates the source array.  */
  index_type sstride[GFC_MAX_DIMENSIONS];
  index_type sstride0;
  index_type srank;
  const GFC_REAL_8 *sptr;

  index_type count[GFC_MAX_DIMENSIONS];
  index_type extent[GFC_MAX_DIMENSIONS];
  index_type n;
  index_type dim;
  index_type ncopies;

  srank = GFC_DESCRIPTOR_RANK(source);

  rrank = srank + 1;
  if (rrank > GFC_MAX_DIMENSIONS)
    runtime_error ("return rank too large in spread()");

  if (along > rrank)
      runtime_error ("dim outside of rank in spread()");

  ncopies = pncopies;

  if (ret->data == NULL)
    {

      size_t ub, stride;

      /* The front end has signalled that we need to populate the
	 return array descriptor.  */
      ret->dtype = (source->dtype & ~GFC_DTYPE_RANK_MASK) | rrank;
      dim = 0;
      rs = 1;
      for (n = 0; n < rrank; n++)
	{
	  stride = rs;
	  if (n == along - 1)
	    {
	      ub = ncopies - 1;
	      rdelta = rs;
	      rs *= ncopies;
	    }
	  else
	    {
	      count[dim] = 0;
	      extent[dim] = GFC_DESCRIPTOR_EXTENT(source,dim);
	      sstride[dim] = GFC_DESCRIPTOR_STRIDE(source,dim);
	      rstride[dim] = rs;

	      ub = extent[dim] - 1;
	      rs *= extent[dim];
	      dim++;
	    }
	  GFC_DIMENSION_SET(ret->dim[n], 0, ub, stride);
	}
      ret->offset = 0;
      if (rs > 0)
        ret->data = internal_malloc_size (rs * sizeof(GFC_REAL_8));
      else
	{
	  ret->data = internal_malloc_size (1);
	  return;
	}
    }
  else
    {
      int zero_sized;

      zero_sized = 0;

      dim = 0;
      if (GFC_DESCRIPTOR_RANK(ret) != rrank)
	runtime_error ("rank mismatch in spread()");

      if (unlikely (compile_options.bounds_check))
	{
	  for (n = 0; n < rrank; n++)
	    {
	      index_type ret_extent;

	      ret_extent = GFC_DESCRIPTOR_EXTENT(ret,n);
	      if (n == along - 1)
		{
		  rdelta = GFC_DESCRIPTOR_STRIDE(ret,n);

		  if (ret_extent != ncopies)
		    runtime_error("Incorrect extent in return value of SPREAD"
				  " intrinsic in dimension %ld: is %ld,"
				  " should be %ld", (long int) n+1,
				  (long int) ret_extent, (long int) ncopies);
		}
	      else
		{
		  count[dim] = 0;
		  extent[dim] = GFC_DESCRIPTOR_EXTENT(source,dim);
		  if (ret_extent != extent[dim])
		    runtime_error("Incorrect extent in return value of SPREAD"
				  " intrinsic in dimension %ld: is %ld,"
				  " should be %ld", (long int) n+1,
				  (long int) ret_extent,
				  (long int) extent[dim]);
		    
		  if (extent[dim] <= 0)
		    zero_sized = 1;
		  sstride[dim] = GFC_DESCRIPTOR_STRIDE(source,dim);
		  rstride[dim] = GFC_DESCRIPTOR_STRIDE(ret,n);
		  dim++;
		}
	    }
	}
      else
	{
	  for (n = 0; n < rrank; n++)
	    {
	      if (n == along - 1)
		{
		  rdelta = GFC_DESCRIPTOR_STRIDE(ret,n);
		}
	      else
		{
		  count[dim] = 0;
		  extent[dim] = GFC_DESCRIPTOR_EXTENT(source,dim);
		  if (extent[dim] <= 0)
		    zero_sized = 1;
		  sstride[dim] = GFC_DESCRIPTOR_STRIDE(source,dim);
		  rstride[dim] = GFC_DESCRIPTOR_STRIDE(ret,n);
		  dim++;
		}
	    }
	}

      if (zero_sized)
	return;

      if (sstride[0] == 0)
	sstride[0] = 1;
    }
  sstride0 = sstride[0];
  rstride0 = rstride[0];
  rptr = ret->data;
  sptr = source->data;

  while (sptr)
    {
      /* Spread this element.  */
      dest = rptr;
      for (n = 0; n < ncopies; n++)
        {
	  *dest = *sptr;
          dest += rdelta;
        }
      /* Advance to the next element.  */
      sptr += sstride0;
      rptr += rstride0;
      count[0]++;
      n = 0;
      while (count[n] == extent[n])
        {
          /* When we get to the end of a dimension, reset it and increment
             the next dimension.  */
          count[n] = 0;
          /* We could precalculate these products, but this is a less
             frequently used path so probably not worth it.  */
          sptr -= sstride[n] * extent[n];
          rptr -= rstride[n] * extent[n];
          n++;
          if (n >= srank)
            {
              /* Break out of the loop.  */
              sptr = NULL;
              break;
            }
          else
            {
              count[n]++;
              sptr += sstride[n];
              rptr += rstride[n];
            }
        }
    }
}

/* This version of spread_internal treats the special case of a scalar
   source.  This is much simpler than the more general case above.  */

void
spread_scalar_r8 (gfc_array_r8 *ret, const GFC_REAL_8 *source,
			const index_type along, const index_type pncopies)
{
  int n;
  int ncopies = pncopies;
  GFC_REAL_8 * restrict dest;
  index_type stride;

  if (GFC_DESCRIPTOR_RANK (ret) != 1)
    runtime_error ("incorrect destination rank in spread()");

  if (along > 1)
    runtime_error ("dim outside of rank in spread()");

  if (ret->data == NULL)
    {
      ret->data = internal_malloc_size (ncopies * sizeof (GFC_REAL_8));
      ret->offset = 0;
      GFC_DIMENSION_SET(ret->dim[0], 0, ncopies - 1, 1);
    }
  else
    {
      if (ncopies - 1 > (GFC_DESCRIPTOR_EXTENT(ret,0) - 1)
			   / GFC_DESCRIPTOR_STRIDE(ret,0))
	runtime_error ("dim too large in spread()");
    }

  dest = ret->data;
  stride = GFC_DESCRIPTOR_STRIDE(ret,0);

  for (n = 0; n < ncopies; n++)
    {
      *dest = *source;
      dest += stride;
    }
}

#endif
Commit	Line	Data
75f2543f	1	/* Special implementation of the SPREAD intrinsic
748086b7	2	Copyright 2008, 2009 Free Software Foundation, Inc.
75f2543f TK	3	Contributed by Thomas Koenig <tkoenig@gcc.gnu.org>, based on
	4	spread_generic.c written by Paul Brook <paul@nowt.org>
	5
	6	This file is part of the GNU Fortran 95 runtime library (libgfortran).
	7
	8	Libgfortran is free software; you can redistribute it and/or
	9	modify it under the terms of the GNU General Public
	10	License as published by the Free Software Foundation; either
748086b7	11	version 3 of the License, or (at your option) any later version.
75f2543f TK	12
	13	Ligbfortran is distributed in the hope that it will be useful,
	14	but WITHOUT ANY WARRANTY; without even the implied warranty of
	15	MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
	16	GNU General Public License for more details.
	17
748086b7 JJ	18	Under Section 7 of GPL version 3, you are granted additional
	19	permissions described in the GCC Runtime Library Exception, version
	20	3.1, as published by the Free Software Foundation.
	21
	22	You should have received a copy of the GNU General Public License and
	23	a copy of the GCC Runtime Library Exception along with this program;
	24	see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
	25	<http://www.gnu.org/licenses/>. */
75f2543f TK	26
	27	#include "libgfortran.h"
	28	#include <stdlib.h>
	29	#include <assert.h>
	30	#include <string.h>
	31
	32
	33	#if defined (HAVE_GFC_REAL_8)
	34
	35	void
	36	spread_r8 (gfc_array_r8 ret, const gfc_array_r8 source,
	37	const index_type along, const index_type pncopies)
	38	{
	39	/* r.* indicates the return array. */
	40	index_type rstride[GFC_MAX_DIMENSIONS];
	41	index_type rstride0;
	42	index_type rdelta = 0;
	43	index_type rrank;
	44	index_type rs;
	45	GFC_REAL_8 *rptr;
5863aacf	46	GFC_REAL_8 * restrict dest;
75f2543f TK	47	/* s.* indicates the source array. */
	48	index_type sstride[GFC_MAX_DIMENSIONS];
	49	index_type sstride0;
	50	index_type srank;
	51	const GFC_REAL_8 *sptr;
	52
	53	index_type count[GFC_MAX_DIMENSIONS];
	54	index_type extent[GFC_MAX_DIMENSIONS];
	55	index_type n;
	56	index_type dim;
	57	index_type ncopies;
	58
	59	srank = GFC_DESCRIPTOR_RANK(source);
	60
	61	rrank = srank + 1;
	62	if (rrank > GFC_MAX_DIMENSIONS)
	63	runtime_error ("return rank too large in spread()");
	64
	65	if (along > rrank)
	66	runtime_error ("dim outside of rank in spread()");
	67
	68	ncopies = pncopies;
	69
	70	if (ret->data == NULL)
	71	{
dfb55fdc TK	72
	73	size_t ub, stride;
	74
75f2543f TK	75	/* The front end has signalled that we need to populate the
	76	return array descriptor. */
	77	ret->dtype = (source->dtype & ~GFC_DTYPE_RANK_MASK) \| rrank;
	78	dim = 0;
	79	rs = 1;
	80	for (n = 0; n < rrank; n++)
	81	{
dfb55fdc	82	stride = rs;
75f2543f TK	83	if (n == along - 1)
75f2543f TK	84	{
dfb55fdc	85	ub = ncopies - 1;
75f2543f TK	86	rdelta = rs;
	87	rs *= ncopies;
	88	}
	89	else
	90	{
	91	count[dim] = 0;
dfb55fdc TK	92	extent[dim] = GFC_DESCRIPTOR_EXTENT(source,dim);
dfb55fdc TK	93	sstride[dim] = GFC_DESCRIPTOR_STRIDE(source,dim);
75f2543f TK	94	rstride[dim] = rs;
75f2543f TK	95
dfb55fdc	96	ub = extent[dim] - 1;
75f2543f TK	97	rs *= extent[dim];
	98	dim++;
	99	}
dfb55fdc	100	GFC_DIMENSION_SET(ret->dim[n], 0, ub, stride);
75f2543f TK	101	}
	102	ret->offset = 0;
	103	if (rs > 0)
	104	ret->data = internal_malloc_size (rs * sizeof(GFC_REAL_8));
	105	else
	106	{
	107	ret->data = internal_malloc_size (1);
	108	return;
	109	}
	110	}
	111	else
	112	{
	113	int zero_sized;
	114
	115	zero_sized = 0;
	116
	117	dim = 0;
	118	if (GFC_DESCRIPTOR_RANK(ret) != rrank)
	119	runtime_error ("rank mismatch in spread()");
	120
9731c4a3	121	if (unlikely (compile_options.bounds_check))
75f2543f TK	122	{
	123	for (n = 0; n < rrank; n++)
	124	{
	125	index_type ret_extent;
	126
dfb55fdc	127	ret_extent = GFC_DESCRIPTOR_EXTENT(ret,n);
75f2543f TK	128	if (n == along - 1)
75f2543f TK	129	{
dfb55fdc	130	rdelta = GFC_DESCRIPTOR_STRIDE(ret,n);
75f2543f TK	131
	132	if (ret_extent != ncopies)
	133	runtime_error("Incorrect extent in return value of SPREAD"
	134	" intrinsic in dimension %ld: is %ld,"
	135	" should be %ld", (long int) n+1,
	136	(long int) ret_extent, (long int) ncopies);
	137	}
	138	else
	139	{
	140	count[dim] = 0;
dfb55fdc	141	extent[dim] = GFC_DESCRIPTOR_EXTENT(source,dim);
75f2543f TK	142	if (ret_extent != extent[dim])
	143	runtime_error("Incorrect extent in return value of SPREAD"
	144	" intrinsic in dimension %ld: is %ld,"
	145	" should be %ld", (long int) n+1,
	146	(long int) ret_extent,
	147	(long int) extent[dim]);
	148
	149	if (extent[dim] <= 0)
	150	zero_sized = 1;
dfb55fdc TK	151	sstride[dim] = GFC_DESCRIPTOR_STRIDE(source,dim);
dfb55fdc TK	152	rstride[dim] = GFC_DESCRIPTOR_STRIDE(ret,n);
75f2543f TK	153	dim++;
	154	}
	155	}
	156	}
	157	else
	158	{
	159	for (n = 0; n < rrank; n++)
	160	{
	161	if (n == along - 1)
	162	{
dfb55fdc	163	rdelta = GFC_DESCRIPTOR_STRIDE(ret,n);
75f2543f TK	164	}
	165	else
	166	{
	167	count[dim] = 0;
dfb55fdc	168	extent[dim] = GFC_DESCRIPTOR_EXTENT(source,dim);
75f2543f TK	169	if (extent[dim] <= 0)
75f2543f TK	170	zero_sized = 1;
dfb55fdc TK	171	sstride[dim] = GFC_DESCRIPTOR_STRIDE(source,dim);
dfb55fdc TK	172	rstride[dim] = GFC_DESCRIPTOR_STRIDE(ret,n);
75f2543f TK	173	dim++;
	174	}
	175	}
	176	}
	177
	178	if (zero_sized)
	179	return;
	180
	181	if (sstride[0] == 0)
	182	sstride[0] = 1;
	183	}
	184	sstride0 = sstride[0];
	185	rstride0 = rstride[0];
	186	rptr = ret->data;
	187	sptr = source->data;
	188
	189	while (sptr)
	190	{
	191	/* Spread this element. */
	192	dest = rptr;
	193	for (n = 0; n < ncopies; n++)
	194	{
	195	dest = sptr;
	196	dest += rdelta;
	197	}
	198	/* Advance to the next element. */
	199	sptr += sstride0;
	200	rptr += rstride0;
	201	count[0]++;
	202	n = 0;
	203	while (count[n] == extent[n])
	204	{
	205	/* When we get to the end of a dimension, reset it and increment
	206	the next dimension. */
	207	count[n] = 0;
	208	/* We could precalculate these products, but this is a less
	209	frequently used path so probably not worth it. */
	210	sptr -= sstride[n] * extent[n];
	211	rptr -= rstride[n] * extent[n];
	212	n++;
	213	if (n >= srank)
	214	{
	215	/* Break out of the loop. */
	216	sptr = NULL;
	217	break;
	218	}
	219	else
	220	{
	221	count[n]++;
	222	sptr += sstride[n];
	223	rptr += rstride[n];
	224	}
	225	}
	226	}
	227	}
	228
	229	/* This version of spread_internal treats the special case of a scalar
	230	source. This is much simpler than the more general case above. */
	231
	232	void
	233	spread_scalar_r8 (gfc_array_r8 ret, const GFC_REAL_8 source,
	234	const index_type along, const index_type pncopies)
	235	{
	236	int n;
237	int ncopies = pncopies;
5863aacf	238	GFC_REAL_8 * restrict dest;
75f2543f TK	239	index_type stride;
	240
	241	if (GFC_DESCRIPTOR_RANK (ret) != 1)
	242	runtime_error ("incorrect destination rank in spread()");
	243
	244	if (along > 1)
	245	runtime_error ("dim outside of rank in spread()");
	246
	247	if (ret->data == NULL)
	248	{
	249	ret->data = internal_malloc_size (ncopies * sizeof (GFC_REAL_8));
	250	ret->offset = 0;
dfb55fdc	251	GFC_DIMENSION_SET(ret->dim[0], 0, ncopies - 1, 1);
75f2543f TK	252	}
	253	else
	254	{
dfb55fdc TK	255	if (ncopies - 1 > (GFC_DESCRIPTOR_EXTENT(ret,0) - 1)
dfb55fdc TK	256	/ GFC_DESCRIPTOR_STRIDE(ret,0))
75f2543f TK	257	runtime_error ("dim too large in spread()");
	258	}
	259
	260	dest = ret->data;
dfb55fdc	261	stride = GFC_DESCRIPTOR_STRIDE(ret,0);
75f2543f TK	262
	263	for (n = 0; n < ncopies; n++)
	264	{
	265	dest = source;
	266	dest += stride;
	267	}
	268	}
	269
	270	#endif
5863aacf	271