old-autovect-branch/libgfortran/generated/maxval_r10.c

   1 /* Implementation of the MAXVAL intrinsic
   2    Copyright 2002 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 2 of the License, or (at your option) any later version.
  11
  12 In addition to the permissions in the GNU General Public License, the
  13 Free Software Foundation gives you unlimited permission to link the
  14 compiled version of this file into combinations with other programs,
  15 and to distribute those combinations without any restriction coming
  16 from the use of this file.  (The General Public License restrictions
  17 do apply in other respects; for example, they cover modification of
  18 the file, and distribution when not linked into a combine
  19 executable.)
  20
  21 Libgfortran is distributed in the hope that it will be useful,
  22 but WITHOUT ANY WARRANTY; without even the implied warranty of
  23 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  24 GNU General Public License for more details.
  25
  26 You should have received a copy of the GNU General Public
  27 License along with libgfortran; see the file COPYING.  If not,
  28 write to the Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
  29 Boston, MA 02110-1301, USA.  */
  30
  31 #include "config.h"
  32 #include <stdlib.h>
  33 #include <assert.h>
  34 #include <float.h>
  35 #include "libgfortran.h"
  36
  37
  38 #if defined (HAVE_GFC_REAL_10) && defined (HAVE_GFC_REAL_10)
  39
  40
  41 extern void maxval_r10 (gfc_array_r10 * const restrict,
  42         gfc_array_r10 * const restrict, const index_type * const restrict);
  43 export_proto(maxval_r10);
  44
  45 void
  46 maxval_r10 (gfc_array_r10 * const restrict retarray,
  47         gfc_array_r10 * const restrict array,
  48         const index_type * const restrict pdim)
  49 {
  50   index_type count[GFC_MAX_DIMENSIONS];
  51   index_type extent[GFC_MAX_DIMENSIONS];
  52   index_type sstride[GFC_MAX_DIMENSIONS];
  53   index_type dstride[GFC_MAX_DIMENSIONS];
  54   const GFC_REAL_10 * restrict base;
  55   GFC_REAL_10 * restrict dest;
  56   index_type rank;
  57   index_type n;
  58   index_type len;
  59   index_type delta;
  60   index_type dim;
  61
  62   /* Make dim zero based to avoid confusion.  */
  63   dim = (*pdim) - 1;
  64   rank = GFC_DESCRIPTOR_RANK (array) - 1;
  65
  66   /* TODO:  It should be a front end job to correctly set the strides.  */
  67
  68   if (array->dim[0].stride == 0)
  69     array->dim[0].stride = 1;
  70
  71   len = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
  72   delta = array->dim[dim].stride;
  73
  74   for (n = 0; n < dim; n++)
  75     {
  76       sstride[n] = array->dim[n].stride;
  77       extent[n] = array->dim[n].ubound + 1 - array->dim[n].lbound;
  78     }
  79   for (n = dim; n < rank; n++)
  80     {
  81       sstride[n] = array->dim[n + 1].stride;
  82       extent[n] =
  83         array->dim[n + 1].ubound + 1 - array->dim[n + 1].lbound;
  84     }
  85
  86   if (retarray->data == NULL)
  87     {
  88       for (n = 0; n < rank; n++)
  89         {
  90           retarray->dim[n].lbound = 0;
  91           retarray->dim[n].ubound = extent[n]-1;
  92           if (n == 0)
  93             retarray->dim[n].stride = 1;
  94           else
  95             retarray->dim[n].stride = retarray->dim[n-1].stride * extent[n-1];
  96         }
  97
  98       retarray->data
  99          = internal_malloc_size (sizeof (GFC_REAL_10)
 100                                  * retarray->dim[rank-1].stride
 101                                  * extent[rank-1]);
 102       retarray->offset = 0;
 103       retarray->dtype = (array->dtype & ~GFC_DTYPE_RANK_MASK) | rank;
 104     }
 105   else
 106     {
 107       if (retarray->dim[0].stride == 0)
 108         retarray->dim[0].stride = 1;
 109
 110       if (rank != GFC_DESCRIPTOR_RANK (retarray))
 111         runtime_error ("rank of return array incorrect");
 112     }
 113
 114   for (n = 0; n < rank; n++)
 115     {
 116       count[n] = 0;
 117       dstride[n] = retarray->dim[n].stride;
 118       if (extent[n] <= 0)
 119         len = 0;
 120     }
 121
 122   base = array->data;
 123   dest = retarray->data;
 124
 125   while (base)
 126     {
 127       const GFC_REAL_10 * restrict src;
 128       GFC_REAL_10 result;
 129       src = base;
 130       {
 131
 132   result = -GFC_REAL_10_HUGE;
 133         if (len <= 0)
 134           *dest = -GFC_REAL_10_HUGE;
 135         else
 136           {
 137             for (n = 0; n < len; n++, src += delta)
 138               {
 139
 140   if (*src > result)
 141     result = *src;
 142           }
 143             *dest = result;
 144           }
 145       }
 146       /* Advance to the next element.  */
 147       count[0]++;
 148       base += sstride[0];
 149       dest += dstride[0];
 150       n = 0;
 151       while (count[n] == extent[n])
 152         {
 153           /* When we get to the end of a dimension, reset it and increment
 154              the next dimension.  */
 155           count[n] = 0;
 156           /* We could precalculate these products, but this is a less
 157              frequently used path so proabably not worth it.  */
 158           base -= sstride[n] * extent[n];
 159           dest -= dstride[n] * extent[n];
 160           n++;
 161           if (n == rank)
 162             {
 163               /* Break out of the look.  */
 164               base = NULL;
 165               break;
 166             }
 167           else
 168             {
 169               count[n]++;
 170               base += sstride[n];
 171               dest += dstride[n];
 172             }
 173         }
 174     }
 175 }
 176
 177
 178 extern void mmaxval_r10 (gfc_array_r10 * const restrict,
 179         gfc_array_r10 * const restrict, const index_type * const restrict,
 180         gfc_array_l4 * const restrict);
 181 export_proto(mmaxval_r10);
 182
 183 void
 184 mmaxval_r10 (gfc_array_r10 * const restrict retarray,
 185         gfc_array_r10 * const restrict array,
 186         const index_type * const restrict pdim,
 187         gfc_array_l4 * const restrict mask)
 188 {
 189   index_type count[GFC_MAX_DIMENSIONS];
 190   index_type extent[GFC_MAX_DIMENSIONS];
 191   index_type sstride[GFC_MAX_DIMENSIONS];
 192   index_type dstride[GFC_MAX_DIMENSIONS];
 193   index_type mstride[GFC_MAX_DIMENSIONS];
 194   GFC_REAL_10 * restrict dest;
 195   const GFC_REAL_10 * restrict base;
 196   const GFC_LOGICAL_4 * restrict mbase;
 197   int rank;
 198   int dim;
 199   index_type n;
 200   index_type len;
 201   index_type delta;
 202   index_type mdelta;
 203
 204   dim = (*pdim) - 1;
 205   rank = GFC_DESCRIPTOR_RANK (array) - 1;
 206
 207   /* TODO:  It should be a front end job to correctly set the strides.  */
 208
 209   if (array->dim[0].stride == 0)
 210     array->dim[0].stride = 1;
 211
 212   if (mask->dim[0].stride == 0)
 213     mask->dim[0].stride = 1;
 214
 215   len = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
 216   if (len <= 0)
 217     return;
 218   delta = array->dim[dim].stride;
 219   mdelta = mask->dim[dim].stride;
 220
 221   for (n = 0; n < dim; n++)
 222     {
 223       sstride[n] = array->dim[n].stride;
 224       mstride[n] = mask->dim[n].stride;
 225       extent[n] = array->dim[n].ubound + 1 - array->dim[n].lbound;
 226     }
 227   for (n = dim; n < rank; n++)
 228     {
 229       sstride[n] = array->dim[n + 1].stride;
 230       mstride[n] = mask->dim[n + 1].stride;
 231       extent[n] =
 232         array->dim[n + 1].ubound + 1 - array->dim[n + 1].lbound;
 233     }
 234
 235   if (retarray->data == NULL)
 236     {
 237       for (n = 0; n < rank; n++)
 238         {
 239           retarray->dim[n].lbound = 0;
 240           retarray->dim[n].ubound = extent[n]-1;
 241           if (n == 0)
 242             retarray->dim[n].stride = 1;
 243           else
 244             retarray->dim[n].stride = retarray->dim[n-1].stride * extent[n-1];
 245         }
 246
 247       retarray->data
 248          = internal_malloc_size (sizeof (GFC_REAL_10)
 249                                  * retarray->dim[rank-1].stride
 250                                  * extent[rank-1]);
 251       retarray->offset = 0;
 252       retarray->dtype = (array->dtype & ~GFC_DTYPE_RANK_MASK) | rank;
 253     }
 254   else
 255     {
 256       if (retarray->dim[0].stride == 0)
 257         retarray->dim[0].stride = 1;
 258
 259       if (rank != GFC_DESCRIPTOR_RANK (retarray))
 260         runtime_error ("rank of return array incorrect");
 261     }
 262
 263   for (n = 0; n < rank; n++)
 264     {
 265       count[n] = 0;
 266       dstride[n] = retarray->dim[n].stride;
 267       if (extent[n] <= 0)
 268         return;
 269     }
 270
 271   dest = retarray->data;
 272   base = array->data;
 273   mbase = mask->data;
 274
 275   if (GFC_DESCRIPTOR_SIZE (mask) != 4)
 276     {
 277       /* This allows the same loop to be used for all logical types.  */
 278       assert (GFC_DESCRIPTOR_SIZE (mask) == 8);
 279       for (n = 0; n < rank; n++)
 280         mstride[n] <<= 1;
 281       mdelta <<= 1;
 282       mbase = (GFOR_POINTER_L8_TO_L4 (mbase));
 283     }
 284
 285   while (base)
 286     {
 287       const GFC_REAL_10 * restrict src;
 288       const GFC_LOGICAL_4 * restrict msrc;
 289       GFC_REAL_10 result;
 290       src = base;
 291       msrc = mbase;
 292       {
 293
 294   result = -GFC_REAL_10_HUGE;
 295         if (len <= 0)
 296           *dest = -GFC_REAL_10_HUGE;
 297         else
 298           {
 299             for (n = 0; n < len; n++, src += delta, msrc += mdelta)
 300               {
 301
 302   if (*msrc && *src > result)
 303     result = *src;
 304               }
 305             *dest = result;
 306           }
 307       }
 308       /* Advance to the next element.  */
 309       count[0]++;
 310       base += sstride[0];
 311       mbase += mstride[0];
 312       dest += dstride[0];
 313       n = 0;
 314       while (count[n] == extent[n])
 315         {
 316           /* When we get to the end of a dimension, reset it and increment
 317              the next dimension.  */
 318           count[n] = 0;
 319           /* We could precalculate these products, but this is a less
 320              frequently used path so proabably not worth it.  */
 321           base -= sstride[n] * extent[n];
 322           mbase -= mstride[n] * extent[n];
 323           dest -= dstride[n] * extent[n];
 324           n++;
 325           if (n == rank)
 326             {
 327               /* Break out of the look.  */
 328               base = NULL;
 329               break;
 330             }
 331           else
 332             {
 333               count[n]++;
 334               base += sstride[n];
 335               mbase += mstride[n];
 336               dest += dstride[n];
 337             }
 338         }
 339     }
 340 }
 341
 342 #endif