libgfortran/generated/sum_r8.c

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