kernel/generic/laswp_ncopy_2.c

   1 /*********************************************************************/
   2 /* Copyright 2009, 2010 The University of Texas at Austin.           */
   3 /* All rights reserved.                                              */
   4 /*                                                                   */
   5 /* Redistribution and use in source and binary forms, with or        */
   6 /* without modification, are permitted provided that the following   */
   7 /* conditions are met:                                               */
   8 /*                                                                   */
   9 /*   1. Redistributions of source code must retain the above         */
  10 /*      copyright notice, this list of conditions and the following  */
  11 /*      disclaimer.                                                  */
  12 /*                                                                   */
  13 /*   2. Redistributions in binary form must reproduce the above      */
  14 /*      copyright notice, this list of conditions and the following  */
  15 /*      disclaimer in the documentation and/or other materials       */
  16 /*      provided with the distribution.                              */
  17 /*                                                                   */
  18 /*    THIS  SOFTWARE IS PROVIDED  BY THE  UNIVERSITY OF  TEXAS AT    */
  19 /*    AUSTIN  ``AS IS''  AND ANY  EXPRESS OR  IMPLIED WARRANTIES,    */
  20 /*    INCLUDING, BUT  NOT LIMITED  TO, THE IMPLIED  WARRANTIES OF    */
  21 /*    MERCHANTABILITY  AND FITNESS FOR  A PARTICULAR  PURPOSE ARE    */
  22 /*    DISCLAIMED.  IN  NO EVENT SHALL THE UNIVERSITY  OF TEXAS AT    */
  23 /*    AUSTIN OR CONTRIBUTORS BE  LIABLE FOR ANY DIRECT, INDIRECT,    */
  24 /*    INCIDENTAL,  SPECIAL, EXEMPLARY,  OR  CONSEQUENTIAL DAMAGES    */
  25 /*    (INCLUDING, BUT  NOT LIMITED TO,  PROCUREMENT OF SUBSTITUTE    */
  26 /*    GOODS  OR  SERVICES; LOSS  OF  USE,  DATA,  OR PROFITS;  OR    */
  27 /*    BUSINESS INTERRUPTION) HOWEVER CAUSED  AND ON ANY THEORY OF    */
  28 /*    LIABILITY, WHETHER  IN CONTRACT, STRICT  LIABILITY, OR TORT    */
  29 /*    (INCLUDING NEGLIGENCE OR OTHERWISE)  ARISING IN ANY WAY OUT    */
  30 /*    OF  THE  USE OF  THIS  SOFTWARE,  EVEN  IF ADVISED  OF  THE    */
  31 /*    POSSIBILITY OF SUCH DAMAGE.                                    */
  32 /*                                                                   */
  33 /* The views and conclusions contained in the software and           */
  34 /* documentation are those of the authors and should not be          */
  35 /* interpreted as representing official policies, either expressed   */
  36 /* or implied, of The University of Texas at Austin.                 */
  37 /*********************************************************************/
  38
  39 #include <stdio.h>
  40 #include "common.h"
  41
  42 #define PREFETCHSIZE 12
  43
  44 #define a2      (a1 + 1)
  45 #define a4      (a3 + 1)
  46
  47 int CNAME(BLASLONG n, BLASLONG k1, BLASLONG k2, FLOAT *a, BLASLONG lda, blasint *ipiv, FLOAT *buffer){
  48
  49   BLASLONG i, j, ip1, ip2;
  50   blasint *piv;
  51   FLOAT *a1, *a3;
  52   FLOAT *b1, *b2, *b3, *b4;
  53   FLOAT A1, A2, B1, B2, A3, A4, B3, B4;
  54
  55   a--;
  56   k1 --;
  57
  58  ipiv += k1;
  59
  60   if (n  <= 0) return 0;
  61
  62   j = (n >> 1);
  63   if (j > 0) {
  64     do {
  65       piv = ipiv;
  66
  67       a1 = a + k1 + 1;
  68       a3 = a1 + 1 * lda;
  69
  70       ip1 = *(piv + 0);
  71       ip2 = *(piv + 1);
  72       piv += 2;
  73
  74       b1 = a + ip1;
  75       b2 = a + ip2;
  76
  77       b3 = b1 + 1 * lda;
  78       b4 = b2 + 1 * lda;
  79
  80       i = ((k2 - k1) >> 1);
  81
  82       if (i > 0) {
  83         do {
  84
  85 #ifdef __GNUC__
  86           __builtin_prefetch(a1 + PREFETCHSIZE, 0, 0);
  87           __builtin_prefetch(a3 + PREFETCHSIZE, 0, 0);
  88 #endif
  89
  90           A1 = *a1;
  91           A2 = *a2;
  92           A3 = *a3;
  93           A4 = *a4;
  94
  95           B1 = *b1;
  96           B2 = *b2;
  97           B3 = *b3;
  98           B4 = *b4;
  99
 100           ip1 = *(piv + 0);
 101           ip2 = *(piv + 1);
 102           piv += 2;
 103
 104           if (b1 == a1) {
 105             if (b2 == a2) {
 106               *(buffer + 0) = A1;
 107               *(buffer + 1) = A3;
 108               *(buffer + 2) = A2;
 109               *(buffer + 3) = A4;
 110             } else {
 111               *(buffer + 0) = A1;
 112               *(buffer + 1) = A3;
 113               *(buffer + 2) = B2;
 114               *(buffer + 3) = B4;
 115
 116               *b2 = A2;
 117               *b4 = A4;
 118             }
 119         } else {
 120           if (b1 == a2) {
 121             if (b2 == a2) {
 122               *(buffer + 0) = A2;
 123               *(buffer + 1) = A4;
 124               *(buffer + 2) = A1;
 125               *(buffer + 3) = A3;
 126             } else {
 127               *(buffer + 0) = A2;
 128               *(buffer + 1) = A4;
 129               *(buffer + 2) = B2;
 130               *(buffer + 3) = B4;
 131               *b2 = A1;
 132               *b4 = A3;
 133             }
 134           } else {
 135             if (b2 == a2) {
 136               *(buffer + 0) = B1;
 137               *(buffer + 1) = B3;
 138               *(buffer + 2) = A2;
 139               *(buffer + 3) = A4;
 140               *b1 = A1;
 141               *b3 = A3;
 142             } else {
 143               if (b2 == b1) {
 144                 *(buffer + 0) = B1;
 145                 *(buffer + 1) = B3;
 146                 *(buffer + 2) = A1;
 147                 *(buffer + 3) = A3;
 148                 *b1 = A2;
 149                 *b3 = A4;
 150               } else {
 151                 *(buffer + 0) = B1;
 152                 *(buffer + 1) = B3;
 153                 *(buffer + 2) = B2;
 154                 *(buffer + 3) = B4;
 155                 *b1 = A1;
 156                 *b2 = A2;
 157                 *b3 = A3;
 158                 *b4 = A4;
 159               }
 160             }
 161           }
 162           }
 163
 164           buffer += 4;
 165
 166           b1 = a + ip1;
 167           b2 = a + ip2;
 168
 169           b3 = b1 + 1 * lda;
 170           b4 = b2 + 1 * lda;
 171
 172           a1 += 2;
 173           a3 += 2;
 174
 175           i --;
 176         } while (i > 0);
 177       }
 178
 179       i = ((k2 - k1) & 1);
 180
 181       if (i > 0) {
 182         A1 = *a1;
 183         B1 = *b1;
 184         A3 = *a3;
 185         B3 = *b3;
 186
 187         if (a1 == b1) {
 188           *(buffer + 0) = A1;
 189           *(buffer + 1) = A3;
 190         } else {
 191           *(buffer + 0) = B1;
 192           *(buffer + 1) = B3;
 193           *b1 = A1;
 194           *b3 = A3;
 195         }
 196         buffer += 2;
 197       }
 198
 199       a += 2 * lda;
 200       j --;
 201     } while (j > 0);
 202   }
 203
 204   if (n & 1) {
 205     piv = ipiv;
 206
 207     a1 = a + k1 + 1;
 208
 209     ip1 = *(piv + 0);
 210     ip2 = *(piv + 1);
 211     piv += 2;
 212
 213     b1 = a + ip1;
 214     b2 = a + ip2;
 215
 216     i = ((k2 - k1) >> 1);
 217
 218     if (i > 0) {
 219       do {
 220         A1 = *a1;
 221         A2 = *a2;
 222         B1 = *b1;
 223         B2 = *b2;
 224
 225         ip1 = *(piv + 0);
 226         ip2 = *(piv + 1);
 227         piv += 2;
 228
 229         if (b1 == a1) {
 230             if (b2 == a2) {
 231               *(buffer + 0) = A1;
 232               *(buffer + 1) = A2;
 233             } else {
 234               *(buffer + 0) = A1;
 235               *(buffer + 1) = B2;
 236
 237               *b2 = A2;
 238             }
 239         } else
 240           if (b1 == a2) {
 241               if (b2 == a2) {
 242                 *(buffer + 0) = A2;
 243                 *(buffer + 1) = A1;
 244               } else {
 245                 *(buffer + 0) = A2;
 246                 *(buffer + 1) = B2;
 247                 *b2 = A1;
 248               }
 249           } else {
 250               if (b2 == a2) {
 251                 *(buffer + 0) = B1;
 252                 *(buffer + 1) = A2;
 253                 *b1 = A1;
 254               } else
 255                 if (b2 == b1) {
 256                   *(buffer + 0) = B1;
 257                   *(buffer + 1) = A1;
 258                   *b1 = A2;
 259                 } else {
 260                   *(buffer + 0) = B1;
 261                   *(buffer + 1) = B2;
 262                   *b1 = A1;
 263                   *b2 = A2;
 264                 }
 265           }
 266
 267         buffer += 2;
 268
 269         b1 = a + ip1;
 270         b2 = a + ip2;
 271
 272         a1 += 2;
 273
 274         i --;
 275       } while (i > 0);
 276     }
 277
 278     i = ((k2 - k1) & 1);
 279
 280     if (i > 0) {
 281       A1 = *a1;
 282       B1 = *b1;
 283
 284       if (a1 == b1) {
 285         *(buffer + 0) = A1;
 286       } else {
 287         *(buffer + 0) = B1;
 288         *b1 = A1;
 289       }
 290     }
 291   }
 292
 293   return 0;
 294 }
 295