libgfortran/generated/matmul_l16.c

   1 /* Implementation of the MATMUL intrinsic
   2    Copyright 2002, 2005, 2006, 2007, 2009 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_LOGICAL_16)
  32
  33 /* Dimensions: retarray(x,y) a(x, count) b(count,y).
  34    Either a or b can be rank 1.  In this case x or y is 1.  */
  35
  36 extern void matmul_l16 (gfc_array_l16 * const restrict,
  37         gfc_array_l1 * const restrict, gfc_array_l1 * const restrict);
  38 export_proto(matmul_l16);
  39
  40 void
  41 matmul_l16 (gfc_array_l16 * const restrict retarray,
  42         gfc_array_l1 * const restrict a, gfc_array_l1 * const restrict b)
  43 {
  44   const GFC_LOGICAL_1 * restrict abase;
  45   const GFC_LOGICAL_1 * restrict bbase;
  46   GFC_LOGICAL_16 * restrict dest;
  47   index_type rxstride;
  48   index_type rystride;
  49   index_type xcount;
  50   index_type ycount;
  51   index_type xstride;
  52   index_type ystride;
  53   index_type x;
  54   index_type y;
  55   int a_kind;
  56   int b_kind;
  57
  58   const GFC_LOGICAL_1 * restrict pa;
  59   const GFC_LOGICAL_1 * restrict pb;
  60   index_type astride;
  61   index_type bstride;
  62   index_type count;
  63   index_type n;
  64
  65   assert (GFC_DESCRIPTOR_RANK (a) == 2
  66           || GFC_DESCRIPTOR_RANK (b) == 2);
  67
  68   if (retarray->data == NULL)
  69     {
  70       if (GFC_DESCRIPTOR_RANK (a) == 1)
  71         {
  72           GFC_DIMENSION_SET(retarray->dim[0], 0,
  73                             GFC_DESCRIPTOR_EXTENT(b,1) - 1, 1);
  74         }
  75       else if (GFC_DESCRIPTOR_RANK (b) == 1)
  76         {
  77           GFC_DIMENSION_SET(retarray->dim[0], 0,
  78                             GFC_DESCRIPTOR_EXTENT(a,0) - 1, 1);
  79         }
  80       else
  81         {
  82           GFC_DIMENSION_SET(retarray->dim[0], 0,
  83                             GFC_DESCRIPTOR_EXTENT(a,0) - 1, 1);
  84
  85           GFC_DIMENSION_SET(retarray->dim[1], 0,
  86                             GFC_DESCRIPTOR_EXTENT(b,1) - 1,
  87                             GFC_DESCRIPTOR_EXTENT(retarray,0));
  88         }
  89
  90       retarray->data
  91         = internal_malloc_size (sizeof (GFC_LOGICAL_16) * size0 ((array_t *) retarray));
  92       retarray->offset = 0;
  93     }
  94     else if (unlikely (compile_options.bounds_check))
  95       {
  96         index_type ret_extent, arg_extent;
  97
  98         if (GFC_DESCRIPTOR_RANK (a) == 1)
  99           {
 100             arg_extent = GFC_DESCRIPTOR_EXTENT(b,1);
 101             ret_extent = GFC_DESCRIPTOR_EXTENT(retarray,0);
 102             if (arg_extent != ret_extent)
 103               runtime_error ("Incorrect extent in return array in"
 104                              " MATMUL intrinsic: is %ld, should be %ld",
 105                              (long int) ret_extent, (long int) arg_extent);
 106           }
 107         else if (GFC_DESCRIPTOR_RANK (b) == 1)
 108           {
 109             arg_extent = GFC_DESCRIPTOR_EXTENT(a,0);
 110             ret_extent = GFC_DESCRIPTOR_EXTENT(retarray,0);
 111             if (arg_extent != ret_extent)
 112               runtime_error ("Incorrect extent in return array in"
 113                              " MATMUL intrinsic: is %ld, should be %ld",
 114                              (long int) ret_extent, (long int) arg_extent);
 115           }
 116         else
 117           {
 118             arg_extent = GFC_DESCRIPTOR_EXTENT(a,0);
 119             ret_extent = GFC_DESCRIPTOR_EXTENT(retarray,0);
 120             if (arg_extent != ret_extent)
 121               runtime_error ("Incorrect extent in return array in"
 122                              " MATMUL intrinsic for dimension 1:"
 123                              " is %ld, should be %ld",
 124                              (long int) ret_extent, (long int) arg_extent);
 125
 126             arg_extent = GFC_DESCRIPTOR_EXTENT(b,1);
 127             ret_extent = GFC_DESCRIPTOR_EXTENT(retarray,1);
 128             if (arg_extent != ret_extent)
 129               runtime_error ("Incorrect extent in return array in"
 130                              " MATMUL intrinsic for dimension 2:"
 131                              " is %ld, should be %ld",
 132                              (long int) ret_extent, (long int) arg_extent);
 133           }
 134       }
 135
 136   abase = a->data;
 137   a_kind = GFC_DESCRIPTOR_SIZE (a);
 138
 139   if (a_kind == 1 || a_kind == 2 || a_kind == 4 || a_kind == 8
 140 #ifdef HAVE_GFC_LOGICAL_16
 141      || a_kind == 16
 142 #endif
 143      )
 144     abase = GFOR_POINTER_TO_L1 (abase, a_kind);
 145   else
 146     internal_error (NULL, "Funny sized logical array");
 147
 148   bbase = b->data;
 149   b_kind = GFC_DESCRIPTOR_SIZE (b);
 150
 151   if (b_kind == 1 || b_kind == 2 || b_kind == 4 || b_kind == 8
 152 #ifdef HAVE_GFC_LOGICAL_16
 153      || b_kind == 16
 154 #endif
 155      )
 156     bbase = GFOR_POINTER_TO_L1 (bbase, b_kind);
 157   else
 158     internal_error (NULL, "Funny sized logical array");
 159
 160   dest = retarray->data;
 161
 162
 163   if (GFC_DESCRIPTOR_RANK (retarray) == 1)
 164     {
 165       rxstride = GFC_DESCRIPTOR_STRIDE(retarray,0);
 166       rystride = rxstride;
 167     }
 168   else
 169     {
 170       rxstride = GFC_DESCRIPTOR_STRIDE(retarray,0);
 171       rystride = GFC_DESCRIPTOR_STRIDE(retarray,1);
 172     }
 173
 174   /* If we have rank 1 parameters, zero the absent stride, and set the size to
 175      one.  */
 176   if (GFC_DESCRIPTOR_RANK (a) == 1)
 177     {
 178       astride = GFC_DESCRIPTOR_STRIDE_BYTES(a,0);
 179       count = GFC_DESCRIPTOR_EXTENT(a,0);
 180       xstride = 0;
 181       rxstride = 0;
 182       xcount = 1;
 183     }
 184   else
 185     {
 186       astride = GFC_DESCRIPTOR_STRIDE_BYTES(a,1);
 187       count = GFC_DESCRIPTOR_EXTENT(a,1);
 188       xstride = GFC_DESCRIPTOR_STRIDE_BYTES(a,0);
 189       xcount = GFC_DESCRIPTOR_EXTENT(a,0);
 190     }
 191   if (GFC_DESCRIPTOR_RANK (b) == 1)
 192     {
 193       bstride = GFC_DESCRIPTOR_STRIDE_BYTES(b,0);
 194       assert(count == GFC_DESCRIPTOR_EXTENT(b,0));
 195       ystride = 0;
 196       rystride = 0;
 197       ycount = 1;
 198     }
 199   else
 200     {
 201       bstride = GFC_DESCRIPTOR_STRIDE_BYTES(b,0);
 202       assert(count == GFC_DESCRIPTOR_EXTENT(b,0));
 203       ystride = GFC_DESCRIPTOR_STRIDE_BYTES(b,1);
 204       ycount = GFC_DESCRIPTOR_EXTENT(b,1);
 205     }
 206
 207   for (y = 0; y < ycount; y++)
 208     {
 209       for (x = 0; x < xcount; x++)
 210         {
 211           /* Do the summation for this element.  For real and integer types
 212              this is the same as DOT_PRODUCT.  For complex types we use do
 213              a*b, not conjg(a)*b.  */
 214           pa = abase;
 215           pb = bbase;
 216           *dest = 0;
 217
 218           for (n = 0; n < count; n++)
 219             {
 220               if (*pa && *pb)
 221                 {
 222                   *dest = 1;
 223                   break;
 224                 }
 225               pa += astride;
 226               pb += bstride;
 227             }
 228
 229           dest += rxstride;
 230           abase += xstride;
 231         }
 232       abase -= xstride * xcount;
 233       bbase += ystride;
 234       dest += rystride - (rxstride * xcount);
 235     }
 236 }
 237
 238 #endif
 239