libgfortran/generated/product_c16.c

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