Optimize daxpy/zaxpy for POWER10
authorRajalakshmi Srinivasaraghavan <rajis@linux.ibm.com>
Thu, 17 Sep 2020 17:56:28 +0000 (12:56 -0500)
committerRajalakshmi Srinivasaraghavan <rajis@linux.ibm.com>
Thu, 17 Sep 2020 17:56:28 +0000 (12:56 -0500)
This patch makes use of new POWER10 vector pair instructions for
loads and stores. Tested in simulator and no new failures.

kernel/power/KERNEL.POWER10
kernel/power/daxpy_microk_power10.c [new file with mode: 0644]
kernel/power/daxpy_power10.c [new file with mode: 0644]
kernel/power/zaxpy_microk_power10.c [new file with mode: 0644]
kernel/power/zaxpy_power10.c [new file with mode: 0644]

index f390fac..ec02e09 100644 (file)
@@ -142,13 +142,13 @@ CASUMKERNEL  = casum.c
 ZASUMKERNEL  = zasum.c
 #
 SAXPYKERNEL  = saxpy.c
-DAXPYKERNEL  = daxpy.c
+DAXPYKERNEL  = daxpy_power10.c
 ifneq ($(GCCVERSIONGTEQ9),1)
 CAXPYKERNEL  = caxpy_power9.S
 else
 CAXPYKERNEL  = caxpy.c
 endif
-ZAXPYKERNEL  = zaxpy.c
+ZAXPYKERNEL  = zaxpy_power10.c
 #
 SCOPYKERNEL  = scopy.c
 DCOPYKERNEL  = dcopy.c
diff --git a/kernel/power/daxpy_microk_power10.c b/kernel/power/daxpy_microk_power10.c
new file mode 100644 (file)
index 0000000..bc9199e
--- /dev/null
@@ -0,0 +1,131 @@
+/***************************************************************************
+Copyright (c) 2020, The OpenBLAS Project
+All rights reserved.
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are
+met:
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+2. Redistributions in binary form must reproduce the above copyright
+notice, this list of conditions and the following disclaimer in
+the documentation and/or other materials provided with the
+distribution.
+3. Neither the name of the OpenBLAS project nor the names of
+its contributors may be used to endorse or promote products
+derived from this software without specific prior written permission.
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+ARE DISCLAIMED. IN NO EVENT SHALL THE OPENBLAS PROJECT OR CONTRIBUTORS BE
+LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
+CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE
+USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+*****************************************************************************/
+
+#define HAVE_KERNEL_8 1
+
+static void daxpy_kernel_8 (long n, double *x, double *y, double alpha)
+{
+  __vector double t0;
+
+  __asm__
+    (
+       XXSPLTD_S(%x4,%x6,0)
+
+       "dcbt           0, %2           \n\t"
+       "dcbt           0, %3           \n\t"
+
+       "lxvp           32, 0(%2)       \n\t"
+       "lxvp           34, 32(%2)      \n\t"
+       "lxvp           40, 64(%2)      \n\t"
+       "lxvp           42, 96(%2)      \n\t"
+
+       "lxvp           36, 0(%3)       \n\t"
+       "lxvp           38, 32(%3)      \n\t"
+       "lxvp           44, 64(%3)      \n\t"
+       "lxvp           46, 96(%3)      \n\t"
+
+       "addi           %2, %2, 128     \n\t"
+
+       "addic.         %1, %1, -16     \n\t"
+       "ble            two%=           \n\t"
+
+       ".align 5                       \n"
+     "one%=:                           \n\t"
+
+       "xvmaddadp      36, 32, %x4     \n\t"
+       "xvmaddadp      37, 33, %x4     \n\t"
+
+       "lxvp           32, 0(%2)       \n\t"
+       "stxvp          36, 0(%3)       \n\t"
+
+       "xvmaddadp      38, 34, %x4     \n\t"
+       "xvmaddadp      39, 35, %x4     \n\t"
+
+       "lxvp           34, 32(%2)      \n\t"
+       "stxvp          38, 32(%3)      \n\t"
+
+
+       "lxvp           36, 128(%3)     \n\t"
+       "lxvp           38, 160(%3)     \n\t"
+
+       "xvmaddadp      44, 40, %x4     \n\t"
+       "xvmaddadp      45, 41, %x4     \n\t"
+
+       "lxvp           40, 64(%2)      \n\t"
+       "stxvp          44, 64(%3)      \n\t"
+
+       "xvmaddadp      46, 42, %x4     \n\t"
+       "xvmaddadp      47, 43, %x4     \n\t"
+
+       "lxvp           42, 96(%2)      \n\t"
+       "stxvp          46, 96(%3)      \n\t"
+
+       "addi           %2, %2, 128     \n\t"
+       "addi           %3, %3, 128     \n\t"
+
+       "lxvp           44, 64(%3)      \n\t"
+       "lxvp           46, 96(%3)      \n\t"
+
+       "addic.         %1, %1, -16     \n\t"
+       "bgt            one%=           \n"
+
+     "two%=:                           \n\t"
+
+       "xvmaddadp      36, 32, %x4     \n\t"
+       "xvmaddadp      37, 33, %x4     \n\t"
+       "xvmaddadp      38, 34, %x4     \n\t"
+       "xvmaddadp      39, 35, %x4     \n\t"
+
+       "xvmaddadp      44, 40, %x4     \n\t"
+       "xvmaddadp      45, 41, %x4     \n\t"
+       "xvmaddadp      46, 42, %x4     \n\t"
+       "xvmaddadp      47, 43, %x4     \n\t"
+
+       "stxvp          36, 0(%3)       \n\t"
+       "stxvp          38, 32(%3)      \n\t"
+       "stxvp          44, 64(%3)      \n\t"
+       "stxvp          46, 96(%3)      \n\t"
+
+     "#n=%1 x=%5=%2 y=%0=%3 alpha=%6 t0=%x4\n"
+     :
+       "+m" (*y),
+       "+r" (n),       // 1
+       "+b" (x),       // 2
+       "+b" (y),       // 3
+       "=wa" (t0)      // 4
+     :
+       "m" (*x),
+       "d" (alpha)     // 6 
+     :
+       "cr0",
+       "vs32","vs33","vs34","vs35","vs36","vs37", "vs38", "vs39",
+       "vs40","vs41","vs42","vs43","vs44","vs45","vs46","vs47"
+     );
+
+}
+
+
diff --git a/kernel/power/daxpy_power10.c b/kernel/power/daxpy_power10.c
new file mode 100644 (file)
index 0000000..ebe91a8
--- /dev/null
@@ -0,0 +1,121 @@
+/***************************************************************************
+Copyright (c) 2020, The OpenBLAS Project
+All rights reserved.
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are
+met:
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+2. Redistributions in binary form must reproduce the above copyright
+notice, this list of conditions and the following disclaimer in
+the documentation and/or other materials provided with the
+distribution.
+3. Neither the name of the OpenBLAS project nor the names of
+its contributors may be used to endorse or promote products
+derived from this software without specific prior written permission.
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+ARE DISCLAIMED. IN NO EVENT SHALL THE OPENBLAS PROJECT OR CONTRIBUTORS BE
+LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
+CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE
+USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+*****************************************************************************/
+
+#include "common.h"
+
+#if defined(__VEC__) || defined(__ALTIVEC__)
+#include "daxpy_microk_power10.c"
+#endif
+
+
+#ifndef HAVE_KERNEL_8
+
+static void daxpy_kernel_8(BLASLONG n, FLOAT *x, FLOAT *y, FLOAT alpha)
+{
+       BLASLONG register i = 0;
+
+       while(i < n)
+        {
+              y[i]   += alpha * x[i];
+              y[i+1] += alpha * x[i+1];
+              y[i+2] += alpha * x[i+2];
+              y[i+3] += alpha * x[i+3];
+              y[i+4] += alpha * x[i+4];
+              y[i+5] += alpha * x[i+5];
+              y[i+6] += alpha * x[i+6];
+              y[i+7] += alpha * x[i+7];
+              i+=8 ;
+
+       }
+
+}
+
+#endif
+
+int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da, FLOAT *x, BLASLONG inc_x, FLOAT *y, BLASLONG inc_y, FLOAT *dummy, BLASLONG dummy2)
+{
+       BLASLONG i=0;
+       BLASLONG ix=0,iy=0;
+
+       if ( n <= 0 )  return(0);
+
+       if ( (inc_x == 1) && (inc_y == 1) )
+       {
+
+               BLASLONG n1 = n & -16;
+
+               if ( n1 )
+                       daxpy_kernel_8(n1, x, y, da);
+
+               i = n1;
+               while(i < n)
+               {
+
+                       y[i] += da * x[i] ;
+                       i++ ;
+
+               }
+               return(0);
+
+
+       }
+
+       BLASLONG n1 = n & -4;
+
+       while(i < n1)
+       {
+
+               FLOAT m1      = da * x[ix] ;
+               FLOAT m2      = da * x[ix+inc_x] ;
+               FLOAT m3      = da * x[ix+2*inc_x] ;
+               FLOAT m4      = da * x[ix+3*inc_x] ;
+
+               y[iy]         += m1 ;
+               y[iy+inc_y]   += m2 ;
+               y[iy+2*inc_y] += m3 ;
+               y[iy+3*inc_y] += m4 ;
+
+               ix  += inc_x*4 ;
+               iy  += inc_y*4 ;
+               i+=4 ;
+
+       }
+
+       while(i < n)
+       {
+
+               y[iy] += da * x[ix] ;
+               ix  += inc_x ;
+               iy  += inc_y ;
+               i++ ;
+
+       }
+       return(0);
+
+}
+
+
diff --git a/kernel/power/zaxpy_microk_power10.c b/kernel/power/zaxpy_microk_power10.c
new file mode 100644 (file)
index 0000000..8e593bb
--- /dev/null
@@ -0,0 +1,200 @@
+/***************************************************************************
+Copyright (c) 2020, The OpenBLAS Project
+All rights reserved.
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are
+met:
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+2. Redistributions in binary form must reproduce the above copyright
+notice, this list of conditions and the following disclaimer in
+the documentation and/or other materials provided with the
+distribution.
+3. Neither the name of the OpenBLAS project nor the names of
+its contributors may be used to endorse or promote products
+derived from this software without specific prior written permission.
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+ARE DISCLAIMED. IN NO EVENT SHALL THE OPENBLAS PROJECT OR CONTRIBUTORS BE
+LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
+CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE
+USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+*****************************************************************************/
+
+#define HAVE_KERNEL_4 1
+static void zaxpy_kernel_4 (long n, double *x, double *y,
+                           double alpha_r, double alpha_i)
+{
+#if !defined(CONJ)
+  static const double mvec[2] = { 1.0, -1.0 };
+#else
+  static const double mvec[2] = { -1.0, 1.0 };
+#endif
+  const double *mvecp = mvec;
+
+  __vector double t0;
+  __vector double t1;
+  __vector double t2;
+  __vector double t3;
+  __vector double t4;
+  __vector double t5;
+  __vector double t6;
+  __vector double t7;
+  long ytmp;
+
+  __asm__
+    (
+       XXSPLTD_S(32,%x15,0)    // alpha_r
+       XXSPLTD_S(33,%x16,0)    // alpha_i
+       "lxvd2x         36, 0, %17      \n\t"   // mvec
+
+#if !defined(CONJ)
+       "xvmuldp                33, 33, 36      \n\t"   // alpha_i * mvec
+#else
+       "xvmuldp                32, 32, 36      \n\t"   // alpha_r * mvec
+#endif
+
+       "mr             %12, %3         \n\t"
+       "dcbt           0, %2           \n\t"
+       "dcbt           0, %3           \n\t"
+
+
+       "lxvp           40, 0(%2)       \n\t"   // x0
+       "lxvp           42, 32(%2)      \n\t"   // x2
+       "lxvp           48, 0(%3)       \n\t"   // y0
+       "lxvp           50, 32(%3)      \n\t"   // y2
+
+       XXSWAPD_S(%x4,40)       // exchange real and imag part
+       XXSWAPD_S(%x5,41)       // exchange real and imag part
+       XXSWAPD_S(%x6,42)       // exchange real and imag part
+       XXSWAPD_S(%x7,43)       // exchange real and imag part
+
+       "lxvp           44, 64(%2)      \n\t"   // x4
+       "lxvp           46, 96(%2)      \n\t"   // x6
+       "lxvp           34, 64(%3)      \n\t"   // y4
+       "lxvp           38, 96(%3)      \n\t"   // y6
+
+       XXSWAPD_S(%x8,44)       // exchange real and imag part
+       XXSWAPD_S(%x9,45)       // exchange real and imag part
+       XXSWAPD_S(%x10,46)      // exchange real and imag part
+       XXSWAPD_S(%x11,47)      // exchange real and imag part
+
+       "addi           %2, %2, 128     \n\t"
+       "addi           %3, %3, 128     \n\t"
+
+       "addic.         %1, %1, -8      \n\t"
+       "ble            two%=           \n\t"
+
+       ".align 5               \n"
+       "one%=:                         \n\t"
+
+       "xvmaddadp      48, 40, 32      \n\t"   // alpha_r * x0_r , alpha_r * x0_i
+       "xvmaddadp      49, 41, 32      \n\t"
+       "lxvp           40, 0(%2)       \n\t"   // x0
+       "xvmaddadp      50, 42, 32      \n\t"
+       "xvmaddadp      51, 43, 32      \n\t"
+       "lxvp           42, 32(%2)      \n\t"   // x2
+
+       "xvmaddadp      34, 44, 32      \n\t"
+       "xvmaddadp      35, 45, 32      \n\t"
+       "lxvp           44, 64(%2)      \n\t"   // x4
+       "xvmaddadp      38, 46, 32      \n\t"
+       "xvmaddadp      39, 47, 32      \n\t"
+       "lxvp           46, 96(%2)      \n\t"   // x6
+
+       "xvmaddadp      48, %x4, 33     \n\t"   // alpha_i * x0_i , alpha_i * x0_r
+       "addi           %2, %2, 128     \n\t"
+       "xvmaddadp      49, %x5, 33     \n\t"
+       "xvmaddadp      50, %x6, 33     \n\t"
+       "xvmaddadp      51, %x7, 33     \n\t"
+
+       "xvmaddadp      34, %x8, 33     \n\t"
+       "xvmaddadp      35, %x9, 33     \n\t"
+       "xvmaddadp      38, %x10, 33    \n\t"
+       "xvmaddadp      39, %x11, 33    \n\t"
+
+       "stxvp          48, 0(%12)      \n\t"
+       "stxvp          50, 32(%12)     \n\t"
+       "stxvp          34, 64(%12)     \n\t"
+       "stxvp          38, 96(%12)     \n\t"
+
+       "addi           %12, %12, 128   \n\t"
+
+       XXSWAPD_S(%x4,40)       // exchange real and imag part
+       XXSWAPD_S(%x5,41)       // exchange real and imag part
+       "lxvp           48, 0(%3)       \n\t"   // y0
+       XXSWAPD_S(%x6,42)       // exchange real and imag part
+       XXSWAPD_S(%x7,43)       // exchange real and imag part
+       "lxvp           50, 32(%3)      \n\t"   // y2
+
+       XXSWAPD_S(%x8,44)       // exchange real and imag part
+       XXSWAPD_S(%x9,45)       // exchange real and imag part
+       "lxvp           34, 64(%3)      \n\t"   // y4
+       XXSWAPD_S(%x10,46)      // exchange real and imag part
+       XXSWAPD_S(%x11,47)      // exchange real and imag part
+       "lxvp           38, 96(%3)      \n\t"   // y6
+
+       "addi           %3, %3, 128     \n\t"
+
+       "addic.         %1, %1, -8      \n\t"
+       "bgt            one%=           \n"
+
+       "two%=:                         \n\t"
+       "xvmaddadp      48, 40, 32      \n\t"   // alpha_r * x0_r , alpha_r * x0_i
+       "xvmaddadp      49, 41, 32      \n\t"
+       "xvmaddadp      50, 42, 32      \n\t"
+       "xvmaddadp      51, 43, 32      \n\t"
+
+       "xvmaddadp      34, 44, 32      \n\t"
+       "xvmaddadp      35, 45, 32      \n\t"
+       "xvmaddadp      38, 46, 32      \n\t"
+       "xvmaddadp      39, 47, 32      \n\t"
+
+       "xvmaddadp      48, %x4, 33     \n\t"   // alpha_i * x0_i , alpha_i * x0_r
+       "xvmaddadp      49, %x5, 33     \n\t"
+       "xvmaddadp      50, %x6, 33     \n\t"
+       "xvmaddadp      51, %x7, 33     \n\t"
+
+       "xvmaddadp      34, %x8, 33     \n\t"
+       "xvmaddadp      35, %x9, 33     \n\t"
+       "xvmaddadp      38, %x10, 33    \n\t"
+       "xvmaddadp      39, %x11, 33    \n\t"
+
+       "stxvp          48, 0(%12)      \n\t"
+       "stxvp          50, 32(%12)     \n\t"
+       "stxvp          34, 64(%12)     \n\t"
+       "stxvp          38, 96(%12)     \n\t"
+
+     "#n=%1 x=%13=%2 y=%0=%3 alpha=(%15,%16) mvecp=%14=%17 ytmp=%12\n"
+     "#t0=%x4 t1=%x5 t2=%x6 t3=%x7 t4=%x8 t5=%x9 t6=%x10 t7=%x11"
+     :
+       "+m" (*y),
+       "+r" (n),       // 1
+       "+b" (x),       // 2
+       "+b" (y),       // 3
+       "=wa" (t0),     // 4
+       "=wa" (t1),     // 5
+       "=wa" (t2),     // 6
+       "=wa" (t3),     // 7
+       "=wa" (t4),     // 8
+       "=wa" (t5),     // 9
+       "=wa" (t6),     // 10
+       "=wa" (t7),     // 11
+       "=b" (ytmp)     // 12
+     :
+       "m" (*x),
+       "m" (*mvecp),
+       "d" (alpha_r),  // 15
+       "d" (alpha_i),  // 16
+       "12" (mvecp)    // 17
+     :
+       "cr0",
+       "vs32","vs33","vs34","vs35","vs36","vs37","vs38","vs39",
+       "vs40","vs41","vs42","vs43","vs44","vs45","vs46","vs47",
+       "vs48","vs49","vs50","vs51"
+     );
+}
diff --git a/kernel/power/zaxpy_power10.c b/kernel/power/zaxpy_power10.c
new file mode 100644 (file)
index 0000000..54cfb8f
--- /dev/null
@@ -0,0 +1,126 @@
+/***************************************************************************
+Copyright (c) 2020, The OpenBLAS Project
+All rights reserved.
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are
+met:
+1. Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+2. Redistributions in binary form must reproduce the above copyright
+notice, this list of conditions and the following disclaimer in
+the documentation and/or other materials provided with the
+distribution.
+3. Neither the name of the OpenBLAS project nor the names of
+its contributors may be used to endorse or promote products
+derived from this software without specific prior written permission.
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+ARE DISCLAIMED. IN NO EVENT SHALL THE OPENBLAS PROJECT OR CONTRIBUTORS BE
+LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
+CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE
+USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+*****************************************************************************/
+
+#include "common.h"
+
+
+#if defined(__VEC__) || defined(__ALTIVEC__)
+#include "zaxpy_microk_power10.c"
+#endif
+
+
+#ifndef HAVE_KERNEL_4
+
+static void zaxpy_kernel_4(BLASLONG n, FLOAT *x, FLOAT *y, FLOAT da_r,FLOAT da_i)
+{
+       BLASLONG register i  = 0;
+       BLASLONG register ix = 0;
+       
+
+       while(i < n)
+        {
+#if !defined(CONJ)
+              y[ix]   += ( da_r * x[ix]   - da_i * x[ix+1] ) ;
+              y[ix+1] += ( da_r * x[ix+1] + da_i * x[ix]   ) ;
+              y[ix+2] += ( da_r * x[ix+2] - da_i * x[ix+3] ) ;
+              y[ix+3] += ( da_r * x[ix+3] + da_i * x[ix+2] ) ;
+#else
+              y[ix]   += ( da_r * x[ix]   + da_i * x[ix+1] ) ;
+              y[ix+1] -= ( da_r * x[ix+1] - da_i * x[ix]   ) ;
+              y[ix+2] += ( da_r * x[ix+2] + da_i * x[ix+3] ) ;
+              y[ix+3] -= ( da_r * x[ix+3] - da_i * x[ix+2] ) ;
+#endif
+
+              ix+=4 ;
+              i+=2 ;
+
+       }
+
+}
+
+#endif
+
+int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, FLOAT *x, BLASLONG inc_x, FLOAT *y, BLASLONG inc_y, FLOAT *dummy, BLASLONG dummy2)
+{
+       BLASLONG i=0;
+       BLASLONG ix=0,iy=0;
+
+       if ( n <= 0 )  return(0);
+
+       if ( (inc_x == 1) && (inc_y == 1) )
+       {
+
+               BLASLONG n1 = n & -16;
+
+               if ( n1 )
+               {
+                       zaxpy_kernel_4 (n1, x, y, da_r, da_i);
+                       ix = 2 * n1;
+               }
+               i = n1;
+               while(i < n)
+               {
+#if !defined(CONJ)
+                       y[ix]   += ( da_r * x[ix]   - da_i * x[ix+1] ) ;
+                       y[ix+1] += ( da_r * x[ix+1] + da_i * x[ix]   ) ;
+#else
+                       y[ix]   += ( da_r * x[ix]   + da_i * x[ix+1] ) ;
+                       y[ix+1] -= ( da_r * x[ix+1] - da_i * x[ix]   ) ;
+#endif
+                       i++ ;
+                       ix += 2;
+
+               }
+               return(0);
+
+
+       }
+
+       inc_x *=2;
+       inc_y *=2;
+
+       while(i < n)
+       {
+
+#if !defined(CONJ)
+                y[iy]   += ( da_r * x[ix]   - da_i * x[ix+1] ) ;
+                y[iy+1] += ( da_r * x[ix+1] + da_i * x[ix]   ) ;
+#else
+                y[iy]   += ( da_r * x[ix]   + da_i * x[ix+1] ) ;
+                y[iy+1] -= ( da_r * x[ix+1] - da_i * x[ix]   ) ;
+#endif
+               ix  += inc_x ;
+               iy  += inc_y ;
+               i++ ;
+
+       }
+       return(0);
+
+}
+
+