libgfortran/m4/spread.m4

   1 `/* Special implementation of the SPREAD intrinsic
   2    Copyright 2008, 2009 Free Software Foundation, Inc.
   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
  11 version 3 of the License, or (at your option) any later version.
  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
  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/>.  */
  26
  27 #include "libgfortran.h"
  28 #include <stdlib.h>
  29 #include <assert.h>
  30 #include <string.h>'
  31
  32 include(iparm.m4)dnl
  33
  34 `#if defined (HAVE_'rtype_name`)
  35
  36 void
  37 spread_'rtype_code` ('rtype` *ret, const 'rtype` *source,
  38                  const index_type along, const index_type pncopies)
  39 {
  40   /* r.* indicates the return array.  */
  41   index_type rstride[GFC_MAX_DIMENSIONS];
  42   index_type rstride0;
  43   index_type rdelta = 0;
  44   index_type rrank;
  45   index_type rs;
  46   'rtype_name` *rptr;
  47   'rtype_name` * restrict dest;
  48   /* s.* indicates the source array.  */
  49   index_type sstride[GFC_MAX_DIMENSIONS];
  50   index_type sstride0;
  51   index_type srank;
  52   const 'rtype_name` *sptr;
  53
  54   index_type count[GFC_MAX_DIMENSIONS];
  55   index_type extent[GFC_MAX_DIMENSIONS];
  56   index_type n;
  57   index_type dim;
  58   index_type ncopies;
  59
  60   srank = GFC_DESCRIPTOR_RANK(source);
  61
  62   rrank = srank + 1;
  63   if (rrank > GFC_MAX_DIMENSIONS)
  64     runtime_error ("return rank too large in spread()");
  65
  66   if (along > rrank)
  67       runtime_error ("dim outside of rank in spread()");
  68
  69   ncopies = pncopies;
  70
  71   if (ret->data == NULL)
  72     {
  73
  74       size_t ub, stride;
  75
  76       /* The front end has signalled that we need to populate the
  77          return array descriptor.  */
  78       ret->dtype = (source->dtype & ~GFC_DTYPE_RANK_MASK) | rrank;
  79       dim = 0;
  80       rs = 1;
  81       for (n = 0; n < rrank; n++)
  82         {
  83           stride = rs;
  84           if (n == along - 1)
  85             {
  86               ub = ncopies - 1;
  87               rdelta = rs;
  88               rs *= ncopies;
  89             }
  90           else
  91             {
  92               count[dim] = 0;
  93               extent[dim] = GFC_DESCRIPTOR_EXTENT(source,dim);
  94               sstride[dim] = GFC_DESCRIPTOR_STRIDE(source,dim);
  95               rstride[dim] = rs;
  96
  97               ub = extent[dim] - 1;
  98               rs *= extent[dim];
  99               dim++;
 100             }
 101           GFC_DIMENSION_SET(ret->dim[n], 0, ub, stride);
 102         }
 103       ret->offset = 0;
 104       if (rs > 0)
 105         ret->data = internal_malloc_size (rs * sizeof('rtype_name`));
 106       else
 107         {
 108           ret->data = internal_malloc_size (1);
 109           return;
 110         }
 111     }
 112   else
 113     {
 114       int zero_sized;
 115
 116       zero_sized = 0;
 117
 118       dim = 0;
 119       if (GFC_DESCRIPTOR_RANK(ret) != rrank)
 120         runtime_error ("rank mismatch in spread()");
 121
 122       if (unlikely (compile_options.bounds_check))
 123         {
 124           for (n = 0; n < rrank; n++)
 125             {
 126               index_type ret_extent;
 127
 128               ret_extent = GFC_DESCRIPTOR_EXTENT(ret,n);
 129               if (n == along - 1)
 130                 {
 131                   rdelta = GFC_DESCRIPTOR_STRIDE(ret,n);
 132
 133                   if (ret_extent != ncopies)
 134                     runtime_error("Incorrect extent in return value of SPREAD"
 135                                   " intrinsic in dimension %ld: is %ld,"
 136                                   " should be %ld", (long int) n+1,
 137                                   (long int) ret_extent, (long int) ncopies);
 138                 }
 139               else
 140                 {
 141                   count[dim] = 0;
 142                   extent[dim] = GFC_DESCRIPTOR_EXTENT(source,dim);
 143                   if (ret_extent != extent[dim])
 144                     runtime_error("Incorrect extent in return value of SPREAD"
 145                                   " intrinsic in dimension %ld: is %ld,"
 146                                   " should be %ld", (long int) n+1,
 147                                   (long int) ret_extent,
 148                                   (long int) extent[dim]);
 149
 150                   if (extent[dim] <= 0)
 151                     zero_sized = 1;
 152                   sstride[dim] = GFC_DESCRIPTOR_STRIDE(source,dim);
 153                   rstride[dim] = GFC_DESCRIPTOR_STRIDE(ret,n);
 154                   dim++;
 155                 }
 156             }
 157         }
 158       else
 159         {
 160           for (n = 0; n < rrank; n++)
 161             {
 162               if (n == along - 1)
 163                 {
 164                   rdelta = GFC_DESCRIPTOR_STRIDE(ret,n);
 165                 }
 166               else
 167                 {
 168                   count[dim] = 0;
 169                   extent[dim] = GFC_DESCRIPTOR_EXTENT(source,dim);
 170                   if (extent[dim] <= 0)
 171                     zero_sized = 1;
 172                   sstride[dim] = GFC_DESCRIPTOR_STRIDE(source,dim);
 173                   rstride[dim] = GFC_DESCRIPTOR_STRIDE(ret,n);
 174                   dim++;
 175                 }
 176             }
 177         }
 178
 179       if (zero_sized)
 180         return;
 181
 182       if (sstride[0] == 0)
 183         sstride[0] = 1;
 184     }
 185   sstride0 = sstride[0];
 186   rstride0 = rstride[0];
 187   rptr = ret->data;
 188   sptr = source->data;
 189
 190   while (sptr)
 191     {
 192       /* Spread this element.  */
 193       dest = rptr;
 194       for (n = 0; n < ncopies; n++)
 195         {
 196           *dest = *sptr;
 197           dest += rdelta;
 198         }
 199       /* Advance to the next element.  */
 200       sptr += sstride0;
 201       rptr += rstride0;
 202       count[0]++;
 203       n = 0;
 204       while (count[n] == extent[n])
 205         {
 206           /* When we get to the end of a dimension, reset it and increment
 207              the next dimension.  */
 208           count[n] = 0;
 209           /* We could precalculate these products, but this is a less
 210              frequently used path so probably not worth it.  */
 211           sptr -= sstride[n] * extent[n];
 212           rptr -= rstride[n] * extent[n];
 213           n++;
 214           if (n >= srank)
 215             {
 216               /* Break out of the loop.  */
 217               sptr = NULL;
 218               break;
 219             }
 220           else
 221             {
 222               count[n]++;
 223               sptr += sstride[n];
 224               rptr += rstride[n];
 225             }
 226         }
 227     }
 228 }
 229
 230 /* This version of spread_internal treats the special case of a scalar
 231    source.  This is much simpler than the more general case above.  */
 232
 233 void
 234 spread_scalar_'rtype_code` ('rtype` *ret, const 'rtype_name` *source,
 235                         const index_type along, const index_type pncopies)
 236 {
 237   int n;
 238   int ncopies = pncopies;
 239   'rtype_name` * restrict dest;
 240   index_type stride;
 241
 242   if (GFC_DESCRIPTOR_RANK (ret) != 1)
 243     runtime_error ("incorrect destination rank in spread()");
 244
 245   if (along > 1)
 246     runtime_error ("dim outside of rank in spread()");
 247
 248   if (ret->data == NULL)
 249     {
 250       ret->data = internal_malloc_size (ncopies * sizeof ('rtype_name`));
 251       ret->offset = 0;
 252       GFC_DIMENSION_SET(ret->dim[0], 0, ncopies - 1, 1);
 253     }
 254   else
 255     {
 256       if (ncopies - 1 > (GFC_DESCRIPTOR_EXTENT(ret,0) - 1)
 257                            / GFC_DESCRIPTOR_STRIDE(ret,0))
 258         runtime_error ("dim too large in spread()");
 259     }
 260
 261   dest = ret->data;
 262   stride = GFC_DESCRIPTOR_STRIDE(ret,0);
 263
 264   for (n = 0; n < ncopies; n++)
 265     {
 266       *dest = *source;
 267       dest += stride;
 268     }
 269 }
 270
 271 #endif
 272 '