libgfortran/generated/product_c16.c

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