libgfortran/generated/maxloc1_4_r8.c

   1 /* Implementation of the MAXLOC 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 (libgfor).
   6
   7 Libgfortran is free software; you can redistribute it and/or
   8 modify it under the terms of the GNU Lesser General Public
   9 License as published by the Free Software Foundation; either
  10 version 2.1 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 Lesser General Public License for more details.
  16
  17 You should have received a copy of the GNU Lesser General Public
  18 License along with libgfor; see the file COPYING.LIB.  If not,
  19 write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
  20 Boston, MA 02111-1307, USA.  */
  21
  22 #include "config.h"
  23 #include <stdlib.h>
  24 #include <assert.h>
  25 #include <float.h>
  26 #include <limits.h>
  27 #include "libgfortran.h"
  28
  29 void
  30 __maxloc1_4_r8 (gfc_array_i4 * retarray, gfc_array_r8 *array, index_type *pdim)
  31 {
  32   index_type count[GFC_MAX_DIMENSIONS - 1];
  33   index_type extent[GFC_MAX_DIMENSIONS - 1];
  34   index_type sstride[GFC_MAX_DIMENSIONS - 1];
  35   index_type dstride[GFC_MAX_DIMENSIONS - 1];
  36   GFC_REAL_8 *base;
  37   GFC_INTEGER_4 *dest;
  38   index_type rank;
  39   index_type n;
  40   index_type len;
  41   index_type delta;
  42   index_type dim;
  43
  44   /* Make dim zero based to avoid confusion.  */
  45   dim = (*pdim) - 1;
  46   rank = GFC_DESCRIPTOR_RANK (array) - 1;
  47   assert (rank == GFC_DESCRIPTOR_RANK (retarray));
  48   if (array->dim[0].stride == 0)
  49     array->dim[0].stride = 1;
  50   if (retarray->dim[0].stride == 0)
  51     retarray->dim[0].stride = 1;
  52
  53   len = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
  54   delta = array->dim[dim].stride;
  55
  56   for (n = 0; n < dim; n++)
  57     {
  58       sstride[n] = array->dim[n].stride;
  59       extent[n] = array->dim[n].ubound + 1 - array->dim[n].lbound;
  60     }
  61   for (n = dim; n < rank; n++)
  62     {
  63       sstride[n] = array->dim[n + 1].stride;
  64       extent[n] =
  65         array->dim[n + 1].ubound + 1 - array->dim[n + 1].lbound;
  66     }
  67
  68   if (retarray->data == NULL)
  69     {
  70       for (n = 0; n < rank; n++)
  71         {
  72           retarray->dim[n].lbound = 0;
  73           retarray->dim[n].ubound = extent[n]-1;
  74           if (n == 0)
  75             retarray->dim[n].stride = 1;
  76           else
  77             retarray->dim[n].stride = retarray->dim[n-1].stride * extent[n-1];
  78         }
  79
  80       retarray->data
  81          = internal_malloc_size (sizeof (GFC_INTEGER_4)
  82                                  * retarray->dim[rank-1].stride
  83                                  * extent[rank-1]);
  84       retarray->base = 0;
  85     }
  86
  87   for (n = 0; n < rank; n++)
  88     {
  89       count[n] = 0;
  90       dstride[n] = retarray->dim[n].stride;
  91       if (extent[n] <= 0)
  92         len = 0;
  93     }
  94
  95   base = array->data;
  96   dest = retarray->data;
  97
  98   while (base)
  99     {
 100       GFC_REAL_8 *src;
 101       GFC_INTEGER_4 result;
 102       src = base;
 103       {
 104
 105   GFC_REAL_8 maxval;
 106   maxval = -GFC_REAL_8_HUGE;
 107   result = 1;
 108         if (len <= 0)
 109           *dest = 0;
 110         else
 111           {
 112             for (n = 0; n < len; n++, src += delta)
 113               {
 114
 115   if (*src > maxval)
 116     {
 117       maxval = *src;
 118       result = (GFC_INTEGER_4)n + 1;
 119     }
 120           }
 121             *dest = result;
 122           }
 123       }
 124       /* Advance to the next element.  */
 125       count[0]++;
 126       base += sstride[0];
 127       dest += dstride[0];
 128       n = 0;
 129       while (count[n] == extent[n])
 130         {
 131           /* When we get to the end of a dimension, reset it and increment
 132              the next dimension.  */
 133           count[n] = 0;
 134           /* We could precalculate these products, but this is a less
 135              frequently used path so proabably not worth it.  */
 136           base -= sstride[n] * extent[n];
 137           dest -= dstride[n] * extent[n];
 138           n++;
 139           if (n == rank)
 140             {
 141               /* Break out of the look.  */
 142               base = NULL;
 143               break;
 144             }
 145           else
 146             {
 147               count[n]++;
 148               base += sstride[n];
 149               dest += dstride[n];
 150             }
 151         }
 152     }
 153 }
 154
 155 void
 156 __mmaxloc1_4_r8 (gfc_array_i4 * retarray, gfc_array_r8 * array, index_type *pdim, gfc_array_l4 * mask)
 157 {
 158   index_type count[GFC_MAX_DIMENSIONS - 1];
 159   index_type extent[GFC_MAX_DIMENSIONS - 1];
 160   index_type sstride[GFC_MAX_DIMENSIONS - 1];
 161   index_type dstride[GFC_MAX_DIMENSIONS - 1];
 162   index_type mstride[GFC_MAX_DIMENSIONS - 1];
 163   GFC_INTEGER_4 *dest;
 164   GFC_REAL_8 *base;
 165   GFC_LOGICAL_4 *mbase;
 166   int rank;
 167   int dim;
 168   index_type n;
 169   index_type len;
 170   index_type delta;
 171   index_type mdelta;
 172
 173   dim = (*pdim) - 1;
 174   rank = GFC_DESCRIPTOR_RANK (array) - 1;
 175   assert (rank == GFC_DESCRIPTOR_RANK (retarray));
 176   if (array->dim[0].stride == 0)
 177     array->dim[0].stride = 1;
 178   if (retarray->dim[0].stride == 0)
 179     retarray->dim[0].stride = 1;
 180
 181   len = array->dim[dim].ubound + 1 - array->dim[dim].lbound;
 182   if (len <= 0)
 183     return;
 184   delta = array->dim[dim].stride;
 185   mdelta = mask->dim[dim].stride;
 186
 187   for (n = 0; n < dim; n++)
 188     {
 189       sstride[n] = array->dim[n].stride;
 190       mstride[n] = mask->dim[n].stride;
 191       extent[n] = array->dim[n].ubound + 1 - array->dim[n].lbound;
 192     }
 193   for (n = dim; n < rank; n++)
 194     {
 195       sstride[n] = array->dim[n + 1].stride;
 196       mstride[n] = mask->dim[n + 1].stride;
 197       extent[n] =
 198         array->dim[n + 1].ubound + 1 - array->dim[n + 1].lbound;
 199     }
 200
 201   for (n = 0; n < rank; n++)
 202     {
 203       count[n] = 0;
 204       dstride[n] = retarray->dim[n].stride;
 205       if (extent[n] <= 0)
 206         return;
 207     }
 208
 209   dest = retarray->data;
 210   base = array->data;
 211   mbase = mask->data;
 212
 213   if (GFC_DESCRIPTOR_SIZE (mask) != 4)
 214     {
 215       /* This allows the same loop to be used for all logical types.  */
 216       assert (GFC_DESCRIPTOR_SIZE (mask) == 8);
 217       for (n = 0; n < rank; n++)
 218         mstride[n] <<= 1;
 219       mdelta <<= 1;
 220       mbase = (GFOR_POINTER_L8_TO_L4 (mbase));
 221     }
 222
 223   while (base)
 224     {
 225       GFC_REAL_8 *src;
 226       GFC_LOGICAL_4 *msrc;
 227       GFC_INTEGER_4 result;
 228       src = base;
 229       msrc = mbase;
 230       {
 231
 232   GFC_REAL_8 maxval;
 233   maxval = -GFC_REAL_8_HUGE;
 234   result = 1;
 235         if (len <= 0)
 236           *dest = 0;
 237         else
 238           {
 239             for (n = 0; n < len; n++, src += delta, msrc += mdelta)
 240               {
 241
 242   if (*msrc && *src > maxval)
 243     {
 244       maxval = *src;
 245       result = (GFC_INTEGER_4)n + 1;
 246     }
 247               }
 248             *dest = result;
 249           }
 250       }
 251       /* Advance to the next element.  */
 252       count[0]++;
 253       base += sstride[0];
 254       mbase += mstride[0];
 255       dest += dstride[0];
 256       n = 0;
 257       while (count[n] == extent[n])
 258         {
 259           /* When we get to the end of a dimension, reset it and increment
 260              the next dimension.  */
 261           count[n] = 0;
 262           /* We could precalculate these products, but this is a less
 263              frequently used path so proabably not worth it.  */
 264           base -= sstride[n] * extent[n];
 265           mbase -= mstride[n] * extent[n];
 266           dest -= dstride[n] * extent[n];
 267           n++;
 268           if (n == rank)
 269             {
 270               /* Break out of the look.  */
 271               base = NULL;
 272               break;
 273             }
 274           else
 275             {
 276               count[n]++;
 277               base += sstride[n];
 278               mbase += mstride[n];
 279               dest += dstride[n];
 280             }
 281         }
 282     }
 283 }
 284