libgfortran/intrinsics/pack_generic.c

   1 /* Generic implementation of the PACK intrinsic
   2    Copyright (C) 2002, 2004, 2005, 2006, 2007 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 Ligbfortran 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 "libgfortran.h"
  32 #include <stdlib.h>
  33 #include <assert.h>
  34 #include <string.h>
  35
  36 /* PACK is specified as follows:
  37
  38    13.14.80 PACK (ARRAY, MASK, [VECTOR])
  39
  40    Description: Pack an array into an array of rank one under the
  41    control of a mask.
  42
  43    Class: Transformational function.
  44
  45    Arguments:
  46       ARRAY   may be of any type. It shall not be scalar.
  47       MASK    shall be of type LOGICAL. It shall be conformable with ARRAY.
  48       VECTOR  (optional) shall be of the same type and type parameters
  49               as ARRAY. VECTOR shall have at least as many elements as
  50               there are true elements in MASK. If MASK is a scalar
  51               with the value true, VECTOR shall have at least as many
  52               elements as there are in ARRAY.
  53
  54    Result Characteristics: The result is an array of rank one with the
  55    same type and type parameters as ARRAY. If VECTOR is present, the
  56    result size is that of VECTOR; otherwise, the result size is the
  57    number /t/ of true elements in MASK unless MASK is scalar with the
  58    value true, in which case the result size is the size of ARRAY.
  59
  60    Result Value: Element /i/ of the result is the element of ARRAY
  61    that corresponds to the /i/th true element of MASK, taking elements
  62    in array element order, for /i/ = 1, 2, ..., /t/. If VECTOR is
  63    present and has size /n/ > /t/, element /i/ of the result has the
  64    value VECTOR(/i/), for /i/ = /t/ + 1, ..., /n/.
  65
  66    Examples: The nonzero elements of an array M with the value
  67    | 0 0 0 |
  68    | 9 0 0 | may be "gathered" by the function PACK. The result of
  69    | 0 0 7 |
  70    PACK (M, MASK = M.NE.0) is [9,7] and the result of PACK (M, M.NE.0,
  71    VECTOR = (/ 2,4,6,8,10,12 /)) is [9,7,6,8,10,12].
  72
  73 There are two variants of the PACK intrinsic: one, where MASK is
  74 array valued, and the other one where MASK is scalar.  */
  75
  76 static void
  77 pack_internal (gfc_array_char *ret, const gfc_array_char *array,
  78                const gfc_array_l1 *mask, const gfc_array_char *vector,
  79                index_type size)
  80 {
  81   /* r.* indicates the return array.  */
  82   index_type rstride0;
  83   char *rptr;
  84   /* s.* indicates the source array.  */
  85   index_type sstride[GFC_MAX_DIMENSIONS];
  86   index_type sstride0;
  87   const char *sptr;
  88   /* m.* indicates the mask array.  */
  89   index_type mstride[GFC_MAX_DIMENSIONS];
  90   index_type mstride0;
  91   const GFC_LOGICAL_1 *mptr;
  92
  93   index_type count[GFC_MAX_DIMENSIONS];
  94   index_type extent[GFC_MAX_DIMENSIONS];
  95   int zero_sized;
  96   index_type n;
  97   index_type dim;
  98   index_type nelem;
  99   index_type total;
 100   int mask_kind;
 101
 102   dim = GFC_DESCRIPTOR_RANK (array);
 103
 104   sptr = array->data;
 105   mptr = mask->data;
 106
 107   /* Use the same loop for all logical types, by using GFC_LOGICAL_1
 108      and using shifting to address size and endian issues.  */
 109
 110   mask_kind = GFC_DESCRIPTOR_SIZE (mask);
 111
 112   if (mask_kind == 1 || mask_kind == 2 || mask_kind == 4 || mask_kind == 8
 113 #ifdef HAVE_GFC_LOGICAL_16
 114       || mask_kind == 16
 115 #endif
 116       )
 117     {
 118       /*  Don't convert a NULL pointer as we use test for NULL below.  */
 119       if (mptr)
 120         mptr = GFOR_POINTER_TO_L1 (mptr, mask_kind);
 121     }
 122   else
 123     runtime_error ("Funny sized logical array");
 124
 125   zero_sized = 0;
 126   for (n = 0; n < dim; n++)
 127     {
 128       count[n] = 0;
 129       extent[n] = array->dim[n].ubound + 1 - array->dim[n].lbound;
 130       if (extent[n] <= 0)
 131        zero_sized = 1;
 132       sstride[n] = array->dim[n].stride * size;
 133       mstride[n] = mask->dim[n].stride * mask_kind;
 134     }
 135   if (sstride[0] == 0)
 136     sstride[0] = size;
 137   if (mstride[0] == 0)
 138     mstride[0] = mask_kind;
 139
 140   if (ret->data == NULL || compile_options.bounds_check)
 141     {
 142       /* Count the elements, either for allocating memory or
 143          for bounds checking.  */
 144
 145       if (vector != NULL)
 146         {
 147           /* The return array will have as many
 148              elements as there are in VECTOR.  */
 149           total = vector->dim[0].ubound + 1 - vector->dim[0].lbound;
 150         }
 151       else
 152         {
 153           /* We have to count the true elements in MASK.  */
 154
 155           /* TODO: We could speed up pack easily in the case of only
 156              few .TRUE. entries in MASK, by keeping track of where we
 157              would be in the source array during the initial traversal
 158              of MASK, and caching the pointers to those elements. Then,
 159              supposed the number of elements is small enough, we would
 160              only have to traverse the list, and copy those elements
 161              into the result array. In the case of datatypes which fit
 162              in one of the integer types we could also cache the
 163              value instead of a pointer to it.
 164              This approach might be bad from the point of view of
 165              cache behavior in the case where our cache is not big
 166              enough to hold all elements that have to be copied.  */
 167
 168           const GFC_LOGICAL_1 *m = mptr;
 169
 170           total = 0;
 171           if (zero_sized)
 172             m = NULL;
 173
 174           while (m)
 175             {
 176               /* Test this element.  */
 177               if (*m)
 178                 total++;
 179
 180               /* Advance to the next element.  */
 181               m += mstride[0];
 182               count[0]++;
 183               n = 0;
 184               while (count[n] == extent[n])
 185                 {
 186                   /* When we get to the end of a dimension, reset it
 187                      and increment the next dimension.  */
 188                   count[n] = 0;
 189                   /* We could precalculate this product, but this is a
 190                      less frequently used path so probably not worth
 191                      it.  */
 192                   m -= mstride[n] * extent[n];
 193                   n++;
 194                   if (n >= dim)
 195                     {
 196                       /* Break out of the loop.  */
 197                       m = NULL;
 198                       break;
 199                     }
 200                   else
 201                     {
 202                       count[n]++;
 203                       m += mstride[n];
 204                     }
 205                 }
 206             }
 207         }
 208
 209       if (ret->data == NULL)
 210         {
 211           /* Setup the array descriptor.  */
 212           ret->dim[0].lbound = 0;
 213           ret->dim[0].ubound = total - 1;
 214           ret->dim[0].stride = 1;
 215
 216           ret->offset = 0;
 217           if (total == 0)
 218             {
 219               /* In this case, nothing remains to be done.  */
 220               ret->data = internal_malloc_size (1);
 221               return;
 222             }
 223           else
 224             ret->data = internal_malloc_size (size * total);
 225         }
 226       else
 227         {
 228           /* We come here because of range checking.  */
 229           index_type ret_extent;
 230
 231           ret_extent = ret->dim[0].ubound + 1 - ret->dim[0].lbound;
 232           if (total != ret_extent)
 233             runtime_error ("Incorrect extent in return value of PACK intrinsic;"
 234                            " is %ld, should be %ld", (long int) total,
 235                            (long int) ret_extent);
 236         }
 237     }
 238
 239   rstride0 = ret->dim[0].stride * size;
 240   if (rstride0 == 0)
 241     rstride0 = size;
 242   sstride0 = sstride[0];
 243   mstride0 = mstride[0];
 244   rptr = ret->data;
 245
 246   while (sptr && mptr)
 247     {
 248       /* Test this element.  */
 249       if (*mptr)
 250         {
 251           /* Add it.  */
 252           memcpy (rptr, sptr, size);
 253           rptr += rstride0;
 254         }
 255       /* Advance to the next element.  */
 256       sptr += sstride0;
 257       mptr += mstride0;
 258       count[0]++;
 259       n = 0;
 260       while (count[n] == extent[n])
 261         {
 262           /* When we get to the end of a dimension, reset it and increment
 263              the next dimension.  */
 264           count[n] = 0;
 265           /* We could precalculate these products, but this is a less
 266              frequently used path so probably not worth it.  */
 267           sptr -= sstride[n] * extent[n];
 268           mptr -= mstride[n] * extent[n];
 269           n++;
 270           if (n >= dim)
 271             {
 272               /* Break out of the loop.  */
 273               sptr = NULL;
 274               break;
 275             }
 276           else
 277             {
 278               count[n]++;
 279               sptr += sstride[n];
 280               mptr += mstride[n];
 281             }
 282         }
 283     }
 284
 285   /* Add any remaining elements from VECTOR.  */
 286   if (vector)
 287     {
 288       n = vector->dim[0].ubound + 1 - vector->dim[0].lbound;
 289       nelem = ((rptr - ret->data) / rstride0);
 290       if (n > nelem)
 291         {
 292           sstride0 = vector->dim[0].stride * size;
 293           if (sstride0 == 0)
 294             sstride0 = size;
 295
 296           sptr = vector->data + sstride0 * nelem;
 297           n -= nelem;
 298           while (n--)
 299             {
 300               memcpy (rptr, sptr, size);
 301               rptr += rstride0;
 302               sptr += sstride0;
 303             }
 304         }
 305     }
 306 }
 307
 308 extern void pack (gfc_array_char *, const gfc_array_char *,
 309                   const gfc_array_l1 *, const gfc_array_char *);
 310 export_proto(pack);
 311
 312 void
 313 pack (gfc_array_char *ret, const gfc_array_char *array,
 314       const gfc_array_l1 *mask, const gfc_array_char *vector)
 315 {
 316   int type;
 317   index_type size;
 318
 319   type = GFC_DESCRIPTOR_TYPE (array);
 320   size = GFC_DESCRIPTOR_SIZE (array);
 321
 322   switch(type)
 323     {
 324     case GFC_DTYPE_INTEGER:
 325     case GFC_DTYPE_LOGICAL:
 326       switch(size)
 327         {
 328         case sizeof (GFC_INTEGER_1):
 329           pack_i1 ((gfc_array_i1 *) ret, (gfc_array_i1 *) array,
 330                    (gfc_array_l1 *) mask, (gfc_array_i1 *) vector);
 331           return;
 332
 333         case sizeof (GFC_INTEGER_2):
 334           pack_i2 ((gfc_array_i2 *) ret, (gfc_array_i2 *) array,
 335                    (gfc_array_l1 *) mask, (gfc_array_i2 *) vector);
 336           return;
 337
 338         case sizeof (GFC_INTEGER_4):
 339           pack_i4 ((gfc_array_i4 *) ret, (gfc_array_i4 *) array,
 340                    (gfc_array_l1 *) mask, (gfc_array_i4 *) vector);
 341           return;
 342
 343         case sizeof (GFC_INTEGER_8):
 344           pack_i8 ((gfc_array_i8 *) ret, (gfc_array_i8 *) array,
 345                    (gfc_array_l1 *) mask, (gfc_array_i8 *) vector);
 346           return;
 347
 348 #ifdef HAVE_GFC_INTEGER_16
 349         case sizeof (GFC_INTEGER_16):
 350           pack_i1 ((gfc_array_i16 *) ret, (gfc_array_i16 *) array,
 351                    (gfc_array_l1 *) mask, (gfc_array_i16 *) vector);
 352           return;
 353 #endif
 354         }
 355     case GFC_DTYPE_REAL:
 356       switch(size)
 357         {
 358         case sizeof (GFC_REAL_4):
 359           pack_r4 ((gfc_array_r4 *) ret, (gfc_array_r4 *) array,
 360                    (gfc_array_l1 *) mask, (gfc_array_r4 *) vector);
 361           return;
 362
 363         case sizeof (GFC_REAL_8):
 364           pack_r8 ((gfc_array_r8 *) ret, (gfc_array_r8 *) array,
 365                    (gfc_array_l1 *) mask, (gfc_array_r8 *) vector);
 366           return;
 367
 368 #ifdef HAVE_GFC_REAL_10
 369         case sizeof (GFC_REAL_10):
 370           pack_r10 ((gfc_array_r10 *) ret, (gfc_array_r10 *) array,
 371                     (gfc_array_l1 *) mask, (gfc_array_r10 *) vector);
 372           return;
 373 #endif
 374
 375 #ifdef HAVE_GFC_REAL_16
 376         case sizeof (GFC_REAL_16):
 377           pack_r16 ((gfc_array_r16 *) ret, (gfc_array_r16 *) array,
 378                     (gfc_array_l1 *) mask, (gfc_array_r16 *) vector);
 379           return;
 380 #endif
 381         }
 382     case GFC_DTYPE_COMPLEX:
 383       switch(size)
 384         {
 385         case sizeof (GFC_COMPLEX_4):
 386           pack_c4 ((gfc_array_c4 *) ret, (gfc_array_c4 *) array,
 387                    (gfc_array_l1 *) mask, (gfc_array_c4 *) vector);
 388           return;
 389
 390         case sizeof (GFC_COMPLEX_8):
 391           pack_c8 ((gfc_array_c8 *) ret, (gfc_array_c8 *) array,
 392                    (gfc_array_l1 *) mask, (gfc_array_c8 *) vector);
 393           return;
 394
 395 #ifdef HAVE_GFC_COMPLEX_10
 396         case sizeof (GFC_COMPLEX_10):
 397           pack_c10 ((gfc_array_c10 *) ret, (gfc_array_c10 *) array,
 398                     (gfc_array_l1 *) mask, (gfc_array_c10 *) vector);
 399           return;
 400 #endif
 401
 402 #ifdef HAVE_GFC_COMPLEX_16
 403         case sizeof (GFC_REAL_16):
 404           pack_c16 ((gfc_array_c16 *) ret, (gfc_array_c16 *) array,
 405                     (gfc_array_l1 *) mask, (gfc_array_c16 *) vector);
 406           return;
 407 #endif
 408
 409         }
 410     }
 411   pack_internal (ret, array, mask, vector, size);
 412 }
 413
 414 extern void pack_char (gfc_array_char *, GFC_INTEGER_4, const gfc_array_char *,
 415                        const gfc_array_l1 *, const gfc_array_char *,
 416                        GFC_INTEGER_4, GFC_INTEGER_4);
 417 export_proto(pack_char);
 418
 419 void
 420 pack_char (gfc_array_char *ret,
 421            GFC_INTEGER_4 ret_length __attribute__((unused)),
 422            const gfc_array_char *array, const gfc_array_l1 *mask,
 423            const gfc_array_char *vector, GFC_INTEGER_4 array_length,
 424            GFC_INTEGER_4 vector_length __attribute__((unused)))
 425 {
 426   pack_internal (ret, array, mask, vector, array_length);
 427 }
 428
 429 static void
 430 pack_s_internal (gfc_array_char *ret, const gfc_array_char *array,
 431                  const GFC_LOGICAL_4 *mask, const gfc_array_char *vector,
 432                  index_type size)
 433 {
 434   /* r.* indicates the return array.  */
 435   index_type rstride0;
 436   char *rptr;
 437   /* s.* indicates the source array.  */
 438   index_type sstride[GFC_MAX_DIMENSIONS];
 439   index_type sstride0;
 440   const char *sptr;
 441
 442   index_type count[GFC_MAX_DIMENSIONS];
 443   index_type extent[GFC_MAX_DIMENSIONS];
 444   index_type n;
 445   index_type dim;
 446   index_type ssize;
 447   index_type nelem;
 448
 449   dim = GFC_DESCRIPTOR_RANK (array);
 450   ssize = 1;
 451   for (n = 0; n < dim; n++)
 452     {
 453       count[n] = 0;
 454       extent[n] = array->dim[n].ubound + 1 - array->dim[n].lbound;
 455       sstride[n] = array->dim[n].stride * size;
 456       ssize *= extent[n];
 457     }
 458   if (sstride[0] == 0)
 459     sstride[0] = size;
 460
 461   sstride0 = sstride[0];
 462   sptr = array->data;
 463
 464   if (ret->data == NULL)
 465     {
 466       /* Allocate the memory for the result.  */
 467       int total;
 468
 469       if (vector != NULL)
 470         {
 471           /* The return array will have as many elements as there are
 472              in vector.  */
 473           total = vector->dim[0].ubound + 1 - vector->dim[0].lbound;
 474         }
 475       else
 476         {
 477           if (*mask)
 478             {
 479               /* The result array will have as many elements as the input
 480                  array.  */
 481               total = extent[0];
 482               for (n = 1; n < dim; n++)
 483                 total *= extent[n];
 484             }
 485           else
 486             /* The result array will be empty.  */
 487             total = 0;
 488         }
 489
 490       /* Setup the array descriptor.  */
 491       ret->dim[0].lbound = 0;
 492       ret->dim[0].ubound = total - 1;
 493       ret->dim[0].stride = 1;
 494       ret->offset = 0;
 495
 496       if (total == 0)
 497         {
 498           ret->data = internal_malloc_size (1);
 499           return;
 500         }
 501       else
 502         ret->data = internal_malloc_size (size * total);
 503     }
 504
 505   rstride0 = ret->dim[0].stride * size;
 506   if (rstride0 == 0)
 507     rstride0 = size;
 508   rptr = ret->data;
 509
 510   /* The remaining possibilities are now:
 511        If MASK is .TRUE., we have to copy the source array into the
 512      result array. We then have to fill it up with elements from VECTOR.
 513        If MASK is .FALSE., we have to copy VECTOR into the result
 514      array. If VECTOR were not present we would have already returned.  */
 515
 516   if (*mask && ssize != 0)
 517     {
 518       while (sptr)
 519         {
 520           /* Add this element.  */
 521           memcpy (rptr, sptr, size);
 522           rptr += rstride0;
 523
 524           /* Advance to the next element.  */
 525           sptr += sstride0;
 526           count[0]++;
 527           n = 0;
 528           while (count[n] == extent[n])
 529             {
 530               /* When we get to the end of a dimension, reset it and
 531                  increment the next dimension.  */
 532               count[n] = 0;
 533               /* We could precalculate these products, but this is a
 534                  less frequently used path so probably not worth it.  */
 535               sptr -= sstride[n] * extent[n];
 536               n++;
 537               if (n >= dim)
 538                 {
 539                   /* Break out of the loop.  */
 540                   sptr = NULL;
 541                   break;
 542                 }
 543               else
 544                 {
 545                   count[n]++;
 546                   sptr += sstride[n];
 547                 }
 548             }
 549         }
 550     }
 551
 552   /* Add any remaining elements from VECTOR.  */
 553   if (vector)
 554     {
 555       n = vector->dim[0].ubound + 1 - vector->dim[0].lbound;
 556       nelem = ((rptr - ret->data) / rstride0);
 557       if (n > nelem)
 558         {
 559           sstride0 = vector->dim[0].stride * size;
 560           if (sstride0 == 0)
 561             sstride0 = size;
 562
 563           sptr = vector->data + sstride0 * nelem;
 564           n -= nelem;
 565           while (n--)
 566             {
 567               memcpy (rptr, sptr, size);
 568               rptr += rstride0;
 569               sptr += sstride0;
 570             }
 571         }
 572     }
 573 }
 574
 575 extern void pack_s (gfc_array_char *ret, const gfc_array_char *array,
 576                     const GFC_LOGICAL_4 *, const gfc_array_char *);
 577 export_proto(pack_s);
 578
 579 void
 580 pack_s (gfc_array_char *ret, const gfc_array_char *array,
 581         const GFC_LOGICAL_4 *mask, const gfc_array_char *vector)
 582 {
 583   pack_s_internal (ret, array, mask, vector, GFC_DESCRIPTOR_SIZE (array));
 584 }
 585
 586 extern void pack_s_char (gfc_array_char *ret, GFC_INTEGER_4,
 587                          const gfc_array_char *array, const GFC_LOGICAL_4 *,
 588                          const gfc_array_char *, GFC_INTEGER_4,
 589                          GFC_INTEGER_4);
 590 export_proto(pack_s_char);
 591
 592 void
 593 pack_s_char (gfc_array_char *ret,
 594              GFC_INTEGER_4 ret_length __attribute__((unused)),
 595              const gfc_array_char *array, const GFC_LOGICAL_4 *mask,
 596              const gfc_array_char *vector, GFC_INTEGER_4 array_length,
 597              GFC_INTEGER_4 vector_length __attribute__((unused)))
 598 {
 599   pack_s_internal (ret, array, mask, vector, array_length);
 600 }