From aa1cebd45b89f48e0d76021f47d3095e8893ba13 Mon Sep 17 00:00:00 2001
From: Chip Kerchner <ckerchner@tenstorrent.com>
Date: Tue, 10 Feb 2026 18:30:02 +0000
Subject: [PATCH] 128-bit versions.

---
 kernel/riscv64/shgemm_kernel_8x8_zvl128b.c | 372 ++++++++++++++++++++-
 1 file changed, 370 insertions(+), 2 deletions(-)

diff --git a/kernel/riscv64/shgemm_kernel_8x8_zvl128b.c b/kernel/riscv64/shgemm_kernel_8x8_zvl128b.c
index 3197f7c00..9e720fd84 100644
--- a/kernel/riscv64/shgemm_kernel_8x8_zvl128b.c
+++ b/kernel/riscv64/shgemm_kernel_8x8_zvl128b.c
@@ -1,4 +1,3 @@
-
 #include "common.h"
 #include <riscv_vector.h>
 
@@ -30,6 +29,16 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat16m1_t A0 = __riscv_vle16_v_f16m1( &A[ai+0*gvl], gvl );
             ai += 8;
 
+#ifdef FP16_NARROW
+            vfloat16m1_t result0 = __riscv_vfmul_vf_f16m1( A0, B0, gvl);
+            vfloat16m1_t result1 = __riscv_vfmul_vf_f16m1( A0, B1, gvl);
+            vfloat16m1_t result2 = __riscv_vfmul_vf_f16m1( A0, B2, gvl);
+            vfloat16m1_t result3 = __riscv_vfmul_vf_f16m1( A0, B3, gvl);
+            vfloat16m1_t result4 = __riscv_vfmul_vf_f16m1( A0, B4, gvl);
+            vfloat16m1_t result5 = __riscv_vfmul_vf_f16m1( A0, B5, gvl);
+            vfloat16m1_t result6 = __riscv_vfmul_vf_f16m1( A0, B6, gvl);
+            vfloat16m1_t result7 = __riscv_vfmul_vf_f16m1( A0, B7, gvl);
+#else
             vfloat32m2_t result0 = __riscv_vfwmul_vf_f32m2( A0, B0, gvl);
             vfloat32m2_t result1 = __riscv_vfwmul_vf_f32m2( A0, B1, gvl);
             vfloat32m2_t result2 = __riscv_vfwmul_vf_f32m2( A0, B2, gvl);
@@ -38,6 +47,7 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat32m2_t result5 = __riscv_vfwmul_vf_f32m2( A0, B5, gvl);
             vfloat32m2_t result6 = __riscv_vfwmul_vf_f32m2( A0, B6, gvl);
             vfloat32m2_t result7 = __riscv_vfwmul_vf_f32m2( A0, B7, gvl);
+#endif
 		
             for(BLASLONG k=1; k<K; k++) {
                 B0 = B[bi+0];
@@ -54,6 +64,16 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 ai += 8;
 
                 
+#ifdef FP16_NARROW
+                result0 = __riscv_vfmacc_vf_f16m1(result0, B0, A0, gvl);
+                result1 = __riscv_vfmacc_vf_f16m1(result1, B1, A0, gvl);
+                result2 = __riscv_vfmacc_vf_f16m1(result2, B2, A0, gvl);
+                result3 = __riscv_vfmacc_vf_f16m1(result3, B3, A0, gvl);
+                result4 = __riscv_vfmacc_vf_f16m1(result4, B4, A0, gvl);
+                result5 = __riscv_vfmacc_vf_f16m1(result5, B5, A0, gvl);
+                result6 = __riscv_vfmacc_vf_f16m1(result6, B6, A0, gvl);
+                result7 = __riscv_vfmacc_vf_f16m1(result7, B7, A0, gvl);
+#else
                 result0 = __riscv_vfwmacc_vf_f32m2(result0, B0, A0, gvl);
                 result1 = __riscv_vfwmacc_vf_f32m2(result1, B1, A0, gvl);
                 result2 = __riscv_vfwmacc_vf_f32m2(result2, B2, A0, gvl);
@@ -62,6 +82,7 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 result5 = __riscv_vfwmacc_vf_f32m2(result5, B5, A0, gvl);
                 result6 = __riscv_vfwmacc_vf_f32m2(result6, B6, A0, gvl);
                 result7 = __riscv_vfwmacc_vf_f32m2(result7, B7, A0, gvl);
+#endif
             }
             
 
@@ -75,7 +96,17 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat32m2_t c5 = __riscv_vle32_v_f32m2( &C[ci], gvl); ci += ldc-gvl*0;
             vfloat32m2_t c6 = __riscv_vle32_v_f32m2( &C[ci], gvl); ci += ldc-gvl*0;
             vfloat32m2_t c7 = __riscv_vle32_v_f32m2( &C[ci], gvl); ci += ldc-gvl*0;
-            
+
+#ifdef FP16_NARROW
+            c0 = __riscv_vfwmacc_vf_f32m2(c0, alpha, result0, gvl);
+            c1 = __riscv_vfwmacc_vf_f32m2(c1, alpha, result1, gvl);
+            c2 = __riscv_vfwmacc_vf_f32m2(c2, alpha, result2, gvl);
+            c3 = __riscv_vfwmacc_vf_f32m2(c3, alpha, result3, gvl);
+            c4 = __riscv_vfwmacc_vf_f32m2(c4, alpha, result4, gvl);
+            c5 = __riscv_vfwmacc_vf_f32m2(c5, alpha, result5, gvl);
+            c6 = __riscv_vfwmacc_vf_f32m2(c6, alpha, result6, gvl);
+            c7 = __riscv_vfwmacc_vf_f32m2(c7, alpha, result7, gvl);
+#else
             c0 = __riscv_vfmacc_vf_f32m2(c0, alpha, result0, gvl);
             c1 = __riscv_vfmacc_vf_f32m2(c1, alpha, result1, gvl);
             c2 = __riscv_vfmacc_vf_f32m2(c2, alpha, result2, gvl);
@@ -84,6 +115,7 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             c5 = __riscv_vfmacc_vf_f32m2(c5, alpha, result5, gvl);
             c6 = __riscv_vfmacc_vf_f32m2(c6, alpha, result6, gvl);
             c7 = __riscv_vfmacc_vf_f32m2(c7, alpha, result7, gvl);
+#endif
 
             ci = n_top * ldc + m_top;
 
@@ -118,6 +150,16 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat16mf2_t A0 = __riscv_vle16_v_f16mf2(&A[ai + 0 * gvl], gvl);
             ai += 4;
 
+#ifdef FP16_NARROW
+            vfloat16mf2_t result0 = __riscv_vfmul_vf_f16mf2( A0, B0, gvl);
+            vfloat16mf2_t result1 = __riscv_vfmul_vf_f16mf2( A0, B1, gvl);
+            vfloat16mf2_t result2 = __riscv_vfmul_vf_f16mf2( A0, B2, gvl);
+            vfloat16mf2_t result3 = __riscv_vfmul_vf_f16mf2( A0, B3, gvl);
+            vfloat16mf2_t result4 = __riscv_vfmul_vf_f16mf2( A0, B4, gvl);
+            vfloat16mf2_t result5 = __riscv_vfmul_vf_f16mf2( A0, B5, gvl);
+            vfloat16mf2_t result6 = __riscv_vfmul_vf_f16mf2( A0, B6, gvl);
+            vfloat16mf2_t result7 = __riscv_vfmul_vf_f16mf2( A0, B7, gvl);
+#else
             vfloat32m1_t result0 = __riscv_vfwmul_vf_f32m1( A0, B0, gvl);
             vfloat32m1_t result1 = __riscv_vfwmul_vf_f32m1( A0, B1, gvl);
             vfloat32m1_t result2 = __riscv_vfwmul_vf_f32m1( A0, B2, gvl);
@@ -126,6 +168,7 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat32m1_t result5 = __riscv_vfwmul_vf_f32m1( A0, B5, gvl);
             vfloat32m1_t result6 = __riscv_vfwmul_vf_f32m1( A0, B6, gvl);
             vfloat32m1_t result7 = __riscv_vfwmul_vf_f32m1( A0, B7, gvl);
+#endif
 
             for(BLASLONG k=1; k < K; ++k) {
                 B0 = B[bi+0];
@@ -141,6 +184,16 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 A0 = __riscv_vle16_v_f16mf2(&A[ai + 0 * gvl], gvl);
                 ai += 4;
 
+#ifdef FP16_NARROW
+                result0 = __riscv_vfmacc_vf_f16mf2(result0, B0, A0, gvl);
+                result1 = __riscv_vfmacc_vf_f16mf2(result1, B1, A0, gvl);
+                result2 = __riscv_vfmacc_vf_f16mf2(result2, B2, A0, gvl);
+                result3 = __riscv_vfmacc_vf_f16mf2(result3, B3, A0, gvl);
+                result4 = __riscv_vfmacc_vf_f16mf2(result4, B4, A0, gvl);
+                result5 = __riscv_vfmacc_vf_f16mf2(result5, B5, A0, gvl);
+                result6 = __riscv_vfmacc_vf_f16mf2(result6, B6, A0, gvl);
+                result7 = __riscv_vfmacc_vf_f16mf2(result7, B7, A0, gvl);
+#else
                 result0 = __riscv_vfwmacc_vf_f32m1(result0, B0, A0, gvl);
                 result1 = __riscv_vfwmacc_vf_f32m1(result1, B1, A0, gvl);
                 result2 = __riscv_vfwmacc_vf_f32m1(result2, B2, A0, gvl);
@@ -149,6 +202,7 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 result5 = __riscv_vfwmacc_vf_f32m1(result5, B5, A0, gvl);
                 result6 = __riscv_vfwmacc_vf_f32m1(result6, B6, A0, gvl);
                 result7 = __riscv_vfwmacc_vf_f32m1(result7, B7, A0, gvl);
+#endif
             }
 
             BLASLONG ci = n_top * ldc + m_top;
@@ -168,6 +222,16 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat32m1_t c6 = __riscv_vle32_v_f32m1(&C[ci], gvl);
             ci += ldc - gvl * 0;
             vfloat32m1_t c7 = __riscv_vle32_v_f32m1(&C[ci], gvl);
+#ifdef FP16_NARROW
+            c0 = __riscv_vfwmacc_vf_f32m1(c0, alpha, result0, gvl);
+            c1 = __riscv_vfwmacc_vf_f32m1(c1, alpha, result1, gvl);
+            c2 = __riscv_vfwmacc_vf_f32m1(c2, alpha, result2, gvl);
+            c3 = __riscv_vfwmacc_vf_f32m1(c3, alpha, result3, gvl);
+            c4 = __riscv_vfwmacc_vf_f32m1(c4, alpha, result4, gvl);
+            c5 = __riscv_vfwmacc_vf_f32m1(c5, alpha, result5, gvl);
+            c6 = __riscv_vfwmacc_vf_f32m1(c6, alpha, result6, gvl);
+            c7 = __riscv_vfwmacc_vf_f32m1(c7, alpha, result7, gvl);
+#else
             c0 = __riscv_vfmacc_vf_f32m1(c0, alpha, result0, gvl);
             c1 = __riscv_vfmacc_vf_f32m1(c1, alpha, result1, gvl);
             c2 = __riscv_vfmacc_vf_f32m1(c2, alpha, result2, gvl);
@@ -176,6 +240,7 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             c5 = __riscv_vfmacc_vf_f32m1(c5, alpha, result5, gvl);
             c6 = __riscv_vfmacc_vf_f32m1(c6, alpha, result6, gvl);
             c7 = __riscv_vfmacc_vf_f32m1(c7, alpha, result7, gvl);
+#endif
 
             ci= n_top * ldc + m_top;
 
@@ -197,6 +262,24 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             BLASLONG ai = m_top * K;
             BLASLONG bi = n_top * K;
 
+#ifdef FP16_NARROW
+            _Float16 result0 = 0;
+            _Float16 result1 = 0;
+            _Float16 result2 = 0;
+            _Float16 result3 = 0;
+            _Float16 result4 = 0;
+            _Float16 result5 = 0;
+            _Float16 result6 = 0;
+            _Float16 result7 = 0;
+            _Float16 result8 = 0;
+            _Float16 result9 = 0;
+            _Float16 result10 = 0;
+            _Float16 result11 = 0;
+            _Float16 result12 = 0;
+            _Float16 result13 = 0;
+            _Float16 result14 = 0;
+            _Float16 result15 = 0;
+#else
             float result0 = 0;
             float result1 = 0;
             float result2 = 0;
@@ -213,8 +296,27 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             float result13 = 0;
             float result14 = 0;
             float result15 = 0;
+#endif
 
             for(BLASLONG k=0; k<K; k++) {
+#ifdef FP16_NARROW
+                result0+=A[ai+0]*B[bi+0];
+                result1+=A[ai+1]*B[bi+0];
+                result2+=A[ai+0]*B[bi+1];
+                result3+=A[ai+1]*B[bi+1];
+                result4+=A[ai+0]*B[bi+2];
+                result5+=A[ai+1]*B[bi+2];
+                result6+=A[ai+0]*B[bi+3];
+                result7+=A[ai+1]*B[bi+3];
+                result8+=A[ai+0]*B[bi+4];
+                result9+=A[ai+1]*B[bi+4];
+                result10+=A[ai+0]*B[bi+5];
+                result11+=A[ai+1]*B[bi+5];
+                result12+=A[ai+0]*B[bi+6];
+                result13+=A[ai+1]*B[bi+6];
+                result14+=A[ai+0]*B[bi+7];
+                result15+=A[ai+1]*B[bi+7];
+#else
                 result0+=(float)(A[ai+0])*(float)(B[bi+0]);
                 result1+=(float)(A[ai+1])*(float)(B[bi+0]);
                 result2+=(float)(A[ai+0])*(float)(B[bi+1]);
@@ -231,12 +333,31 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 result13+=(float)(A[ai+1])*(float)(B[bi+6]);
                 result14+=(float)(A[ai+0])*(float)(B[bi+7]);
                 result15+=(float)(A[ai+1])*(float)(B[bi+7]);
+#endif
                 ai+=2;
                 bi+=8;
             }
             
             
             BLASLONG ci=n_top*ldc+m_top;
+#ifdef FP16_NARROW
+            C[ci + 0 * ldc + 0] += alpha * (float)(result0);
+            C[ci + 0 * ldc + 1] += alpha * (float)(result1);
+            C[ci + 1 * ldc + 0] += alpha * (float)(result2);
+            C[ci + 1 * ldc + 1] += alpha * (float)(result3);
+            C[ci + 2 * ldc + 0] += alpha * (float)(result4);
+            C[ci + 2 * ldc + 1] += alpha * (float)(result5);
+            C[ci + 3 * ldc + 0] += alpha * (float)(result6);
+            C[ci + 3 * ldc + 1] += alpha * (float)(result7);
+            C[ci + 4 * ldc + 0] += alpha * (float)(result8);
+            C[ci + 4 * ldc + 1] += alpha * (float)(result9);
+            C[ci + 5 * ldc + 0] += alpha * (float)(result10);
+            C[ci + 5 * ldc + 1] += alpha * (float)(result11);
+            C[ci + 6 * ldc + 0] += alpha * (float)(result12);
+            C[ci + 6 * ldc + 1] += alpha * (float)(result13);
+            C[ci + 7 * ldc + 0] += alpha * (float)(result14);
+            C[ci + 7 * ldc + 1] += alpha * (float)(result15);
+#else
             C[ci + 0 * ldc + 0] += alpha * result0;
             C[ci + 0 * ldc + 1] += alpha * result1;
             C[ci + 1 * ldc + 0] += alpha * result2;
@@ -253,6 +374,7 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             C[ci + 6 * ldc + 1] += alpha * result13;
             C[ci + 7 * ldc + 0] += alpha * result14;
             C[ci + 7 * ldc + 1] += alpha * result15;
+#endif
 
             m_top+=2;
         }
@@ -260,6 +382,16 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
 
         if( M & 1 ) {
             
+#ifdef FP16_NARROW
+            _Float16 result0 = 0;
+            _Float16 result1 = 0;
+            _Float16 result2 = 0;
+            _Float16 result3 = 0;
+            _Float16 result4 = 0;
+            _Float16 result5 = 0;
+            _Float16 result6 = 0;
+            _Float16 result7 = 0;
+#else
             float result0 = 0;
             float result1 = 0;
             float result2 = 0;
@@ -268,11 +400,22 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             float result5 = 0;
             float result6 = 0;
             float result7 = 0;
+#endif
             
             BLASLONG ai = m_top * K;
             BLASLONG bi = n_top * K;
 
             for(BLASLONG k=0; k<K; k++) {
+#ifdef FP16_NARROW
+                result0+=A[ai+0]*B[bi+0];
+                result1+=A[ai+0]*B[bi+1];
+                result2+=A[ai+0]*B[bi+2];
+                result3+=A[ai+0]*B[bi+3];
+                result4+=A[ai+0]*B[bi+4];
+                result5+=A[ai+0]*B[bi+5];
+                result6+=A[ai+0]*B[bi+6];
+                result7+=A[ai+0]*B[bi+7];
+#else
                 result0+=(float)(A[ai+0])*(float)(B[bi+0]);
                 result1+=(float)(A[ai+0])*(float)(B[bi+1]);
                 result2+=(float)(A[ai+0])*(float)(B[bi+2]);
@@ -281,11 +424,22 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 result5+=(float)(A[ai+0])*(float)(B[bi+5]);
                 result6+=(float)(A[ai+0])*(float)(B[bi+6]);
                 result7+=(float)(A[ai+0])*(float)(B[bi+7]);
+#endif
                 ai+=1;
                 bi+=8;
             }
 
             BLASLONG ci = n_top * ldc + m_top;
+#ifdef FP16_NARROW
+            C[ci + 0 * ldc + 0] += alpha * (float)(result0);
+            C[ci + 1 * ldc + 0] += alpha * (float)(result1);
+            C[ci + 2 * ldc + 0] += alpha * (float)(result2);
+            C[ci + 3 * ldc + 0] += alpha * (float)(result3);
+            C[ci + 4 * ldc + 0] += alpha * (float)(result4);
+            C[ci + 5 * ldc + 0] += alpha * (float)(result5);
+            C[ci + 6 * ldc + 0] += alpha * (float)(result6);
+            C[ci + 7 * ldc + 0] += alpha * (float)(result7);
+#else
             C[ci + 0 * ldc + 0] += alpha * result0;
             C[ci + 1 * ldc + 0] += alpha * result1;
             C[ci + 2 * ldc + 0] += alpha * result2;
@@ -294,6 +448,7 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             C[ci + 5 * ldc + 0] += alpha * result5;
             C[ci + 6 * ldc + 0] += alpha * result6;
             C[ci + 7 * ldc + 0] += alpha * result7;
+#endif
             m_top+=1;
         }
 
@@ -318,10 +473,17 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat16m1_t A0 = __riscv_vle16_v_f16m1( &A[ai+0*gvl], gvl );
             ai += 8;
 
+#ifdef FP16_NARROW
+            vfloat16m1_t result0 = __riscv_vfmul_vf_f16m1( A0, B0, gvl);
+            vfloat16m1_t result1 = __riscv_vfmul_vf_f16m1( A0, B1, gvl);
+            vfloat16m1_t result2 = __riscv_vfmul_vf_f16m1( A0, B2, gvl);
+            vfloat16m1_t result3 = __riscv_vfmul_vf_f16m1( A0, B3, gvl);
+#else
             vfloat32m2_t result0 = __riscv_vfwmul_vf_f32m2( A0, B0, gvl);
             vfloat32m2_t result1 = __riscv_vfwmul_vf_f32m2( A0, B1, gvl);
             vfloat32m2_t result2 = __riscv_vfwmul_vf_f32m2( A0, B2, gvl);
             vfloat32m2_t result3 = __riscv_vfwmul_vf_f32m2( A0, B3, gvl);
+#endif
 		
             for(BLASLONG k=1; k<K; k++) {
                 B0 = B[bi+0];
@@ -333,10 +495,17 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 A0 = __riscv_vle16_v_f16m1( &A[ai+0*gvl], gvl );
                 ai += 8;
 
+#ifdef FP16_NARROW
+                result0 = __riscv_vfmacc_vf_f16m1(result0, B0, A0, gvl);
+                result1 = __riscv_vfmacc_vf_f16m1(result1, B1, A0, gvl);
+                result2 = __riscv_vfmacc_vf_f16m1(result2, B2, A0, gvl);
+                result3 = __riscv_vfmacc_vf_f16m1(result3, B3, A0, gvl);
+#else
                 result0 = __riscv_vfwmacc_vf_f32m2(result0, B0, A0, gvl);
                 result1 = __riscv_vfwmacc_vf_f32m2(result1, B1, A0, gvl);
                 result2 = __riscv_vfwmacc_vf_f32m2(result2, B2, A0, gvl);
                 result3 = __riscv_vfwmacc_vf_f32m2(result3, B3, A0, gvl);
+#endif
             }
             
 
@@ -347,10 +516,17 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat32m2_t c2 = __riscv_vle32_v_f32m2( &C[ci], gvl); ci += ldc - gvl * 0;
             vfloat32m2_t c3 = __riscv_vle32_v_f32m2( &C[ci], gvl);
             
+#ifdef FP16_NARROW
+            c0 = __riscv_vfwmacc_vf_f32m2(c0, alpha, result0, gvl);
+            c1 = __riscv_vfwmacc_vf_f32m2(c1, alpha, result1, gvl);
+            c2 = __riscv_vfwmacc_vf_f32m2(c2, alpha, result2, gvl);
+            c3 = __riscv_vfwmacc_vf_f32m2(c3, alpha, result3, gvl);
+#else
             c0 = __riscv_vfmacc_vf_f32m2(c0, alpha, result0, gvl);
             c1 = __riscv_vfmacc_vf_f32m2(c1, alpha, result1, gvl);
             c2 = __riscv_vfmacc_vf_f32m2(c2, alpha, result2, gvl);
             c3 = __riscv_vfmacc_vf_f32m2(c3, alpha, result3, gvl);
+#endif
 
             ci = n_top * ldc + m_top;
 
@@ -375,10 +551,17 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat16mf2_t A0 = __riscv_vle16_v_f16mf2(&A[ai + 0 * gvl], gvl);
             ai += 4;
 
+#ifdef FP16_NARROW
+            vfloat16mf2_t result0 = __riscv_vfmul_vf_f16mf2( A0, B0, gvl);
+            vfloat16mf2_t result1 = __riscv_vfmul_vf_f16mf2( A0, B1, gvl);
+            vfloat16mf2_t result2 = __riscv_vfmul_vf_f16mf2( A0, B2, gvl);
+            vfloat16mf2_t result3 = __riscv_vfmul_vf_f16mf2( A0, B3, gvl);
+#else
             vfloat32m1_t result0 = __riscv_vfwmul_vf_f32m1( A0, B0, gvl);
             vfloat32m1_t result1 = __riscv_vfwmul_vf_f32m1( A0, B1, gvl);
             vfloat32m1_t result2 = __riscv_vfwmul_vf_f32m1( A0, B2, gvl);
             vfloat32m1_t result3 = __riscv_vfwmul_vf_f32m1( A0, B3, gvl);
+#endif
 
             for(BLASLONG k=1; k < K; ++k) {
                 B0 = B[bi+0];
@@ -390,10 +573,17 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 A0 = __riscv_vle16_v_f16mf2(&A[ai + 0 * gvl], gvl);
                 ai += 4;
 
+#ifdef FP16_NARROW
+                result0 = __riscv_vfmacc_vf_f16mf2(result0, B0, A0, gvl);
+                result1 = __riscv_vfmacc_vf_f16mf2(result1, B1, A0, gvl);
+                result2 = __riscv_vfmacc_vf_f16mf2(result2, B2, A0, gvl);
+                result3 = __riscv_vfmacc_vf_f16mf2(result3, B3, A0, gvl);
+#else
                 result0 = __riscv_vfwmacc_vf_f32m1(result0, B0, A0, gvl);
                 result1 = __riscv_vfwmacc_vf_f32m1(result1, B1, A0, gvl);
                 result2 = __riscv_vfwmacc_vf_f32m1(result2, B2, A0, gvl);
                 result3 = __riscv_vfwmacc_vf_f32m1(result3, B3, A0, gvl);
+#endif
             }
 
             BLASLONG ci = n_top * ldc + m_top;
@@ -405,10 +595,17 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat32m1_t c2 = __riscv_vle32_v_f32m1(&C[ci], gvl);
             ci += ldc - gvl * 0;
             vfloat32m1_t c3 = __riscv_vle32_v_f32m1(&C[ci], gvl);
+#ifdef FP16_NARROW
+            c0 = __riscv_vfwmacc_vf_f32m1(c0, alpha, result0, gvl);
+            c1 = __riscv_vfwmacc_vf_f32m1(c1, alpha, result1, gvl);
+            c2 = __riscv_vfwmacc_vf_f32m1(c2, alpha, result2, gvl);
+            c3 = __riscv_vfwmacc_vf_f32m1(c3, alpha, result3, gvl);
+#else
             c0 = __riscv_vfmacc_vf_f32m1(c0, alpha, result0, gvl);
             c1 = __riscv_vfmacc_vf_f32m1(c1, alpha, result1, gvl);
             c2 = __riscv_vfmacc_vf_f32m1(c2, alpha, result2, gvl);
             c3 = __riscv_vfmacc_vf_f32m1(c3, alpha, result3, gvl);
+#endif
 
             ci= n_top * ldc + m_top;
 
@@ -426,6 +623,16 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             BLASLONG ai = m_top * K;
             BLASLONG bi = n_top * K;
 
+#ifdef FP16_NARROW
+            _Float16 result0 = 0;
+            _Float16 result1 = 0;
+            _Float16 result2 = 0;
+            _Float16 result3 = 0;
+            _Float16 result4 = 0;
+            _Float16 result5 = 0;
+            _Float16 result6 = 0;
+            _Float16 result7 = 0;
+#else
             float result0 = 0;
             float result1 = 0;
             float result2 = 0;
@@ -434,8 +641,19 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             float result5 = 0;
             float result6 = 0;
             float result7 = 0;
+#endif
 
             for(BLASLONG k=0; k<K; k++) {
+#ifdef FP16_NARROW
+                result0+=A[ai+0]*B[bi+0];
+                result1+=A[ai+1]*B[bi+0];
+                result2+=A[ai+0]*B[bi+1];
+                result3+=A[ai+1]*B[bi+1];
+                result4+=A[ai+0]*B[bi+2];
+                result5+=A[ai+1]*B[bi+2];
+                result6+=A[ai+0]*B[bi+3];
+                result7+=A[ai+1]*B[bi+3];
+#else
                 result0+=(float)(A[ai+0])*(float)(B[bi+0]);
                 result1+=(float)(A[ai+1])*(float)(B[bi+0]);
                 result2+=(float)(A[ai+0])*(float)(B[bi+1]);
@@ -444,12 +662,23 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 result5+=(float)(A[ai+1])*(float)(B[bi+2]);
                 result6+=(float)(A[ai+0])*(float)(B[bi+3]);
                 result7+=(float)(A[ai+1])*(float)(B[bi+3]);
+#endif
                 ai+=2;
                 bi+=4;
             }
             
             
             BLASLONG ci=n_top*ldc+m_top;
+#ifdef FP16_NARROW
+            C[ci + 0 * ldc + 0] += alpha * (float)(result0);
+            C[ci + 0 * ldc + 1] += alpha * (float)(result1);
+            C[ci + 1 * ldc + 0] += alpha * (float)(result2);
+            C[ci + 1 * ldc + 1] += alpha * (float)(result3);
+            C[ci + 2 * ldc + 0] += alpha * (float)(result4);
+            C[ci + 2 * ldc + 1] += alpha * (float)(result5);
+            C[ci + 3 * ldc + 0] += alpha * (float)(result6);
+            C[ci + 3 * ldc + 1] += alpha * (float)(result7);
+#else
             C[ci + 0 * ldc + 0] += alpha * result0;
             C[ci + 0 * ldc + 1] += alpha * result1;
             C[ci + 1 * ldc + 0] += alpha * result2;
@@ -458,6 +687,7 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             C[ci + 2 * ldc + 1] += alpha * result5;
             C[ci + 3 * ldc + 0] += alpha * result6;
             C[ci + 3 * ldc + 1] += alpha * result7;
+#endif
 
             m_top += 2;
         }
@@ -465,28 +695,49 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
 
         if( M & 1 ) {
             
+#ifdef FP16_NARROW
+            _Float16 result0 = 0;
+            _Float16 result1 = 0;
+            _Float16 result2 = 0;
+            _Float16 result3 = 0;
+#else
             float result0 = 0;
             float result1 = 0;
             float result2 = 0;
             float result3 = 0;
+#endif
             
             BLASLONG ai = m_top * K;
             BLASLONG bi = n_top * K;
 
             for(BLASLONG k=0; k<K; k++) {
+#ifdef FP16_NARROW
+                result0+=A[ai+0]*B[bi+0];
+                result1+=A[ai+0]*B[bi+1];
+                result2+=A[ai+0]*B[bi+2];
+                result3+=A[ai+0]*B[bi+3];
+#else
                 result0+=(float)(A[ai+0])*(float)(B[bi+0]);
                 result1+=(float)(A[ai+0])*(float)(B[bi+1]);
                 result2+=(float)(A[ai+0])*(float)(B[bi+2]);
                 result3+=(float)(A[ai+0])*(float)(B[bi+3]);
+#endif
                 ai+=1;
                 bi+=4;
             }
 
             BLASLONG ci = n_top * ldc + m_top;
+#ifdef FP16_NARROW
+            C[ci + 0 * ldc + 0] += alpha * (float)(result0);
+            C[ci + 1 * ldc + 0] += alpha * (float)(result1);
+            C[ci + 2 * ldc + 0] += alpha * (float)(result2);
+            C[ci + 3 * ldc + 0] += alpha * (float)(result3);
+#else
             C[ci + 0 * ldc + 0] += alpha * result0;
             C[ci + 1 * ldc + 0] += alpha * result1;
             C[ci + 2 * ldc + 0] += alpha * result2;
             C[ci + 3 * ldc + 0] += alpha * result3;
+#endif
             m_top += 1;
         }
 
@@ -511,8 +762,13 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat16m1_t A0 = __riscv_vle16_v_f16m1( &A[ai+0*gvl], gvl );
             ai += 8;
 
+#ifdef FP16_NARROW
+            vfloat16m1_t result0 = __riscv_vfmul_vf_f16m1( A0, B0, gvl);
+            vfloat16m1_t result1 = __riscv_vfmul_vf_f16m1( A0, B1, gvl);
+#else
             vfloat32m2_t result0 = __riscv_vfwmul_vf_f32m2( A0, B0, gvl);
             vfloat32m2_t result1 = __riscv_vfwmul_vf_f32m2( A0, B1, gvl);
+#endif
 		
             for(BLASLONG k=1; k<K; k++) {
                 B0 = B[bi+0];
@@ -522,8 +778,13 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 A0 = __riscv_vle16_v_f16m1( &A[ai+0*gvl], gvl );
                 ai += 8;
 
+#ifdef FP16_NARROW
+                result0 = __riscv_vfmacc_vf_f16m1(result0, B0, A0, gvl);
+                result1 = __riscv_vfmacc_vf_f16m1(result1, B1, A0, gvl);
+#else
                 result0 = __riscv_vfwmacc_vf_f32m2(result0, B0, A0, gvl);
                 result1 = __riscv_vfwmacc_vf_f32m2(result1, B1, A0, gvl);
+#endif
             }
 
 
@@ -532,8 +793,13 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat32m2_t c0 = __riscv_vle32_v_f32m2( &C[ci], gvl); ci += ldc - gvl * 0;
             vfloat32m2_t c1 = __riscv_vle32_v_f32m2( &C[ci], gvl); 
             
+#ifdef FP16_NARROW
+            c0 = __riscv_vfwmacc_vf_f32m2(c0, alpha, result0, gvl);
+            c1 = __riscv_vfwmacc_vf_f32m2(c1, alpha, result1, gvl);
+#else
             c0 = __riscv_vfmacc_vf_f32m2(c0, alpha, result0, gvl);
             c1 = __riscv_vfmacc_vf_f32m2(c1, alpha, result1, gvl);
+#endif
 
             ci = n_top * ldc + m_top;
 
@@ -554,8 +820,13 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat16mf2_t A0 = __riscv_vle16_v_f16mf2(&A[ai + 0 * gvl], gvl);
             ai += 4;
 
+#ifdef FP16_NARROW
+            vfloat16mf2_t result0 = __riscv_vfmul_vf_f16mf2( A0, B0, gvl);
+            vfloat16mf2_t result1 = __riscv_vfmul_vf_f16mf2( A0, B1, gvl);
+#else
             vfloat32m1_t result0 = __riscv_vfwmul_vf_f32m1( A0, B0, gvl);
             vfloat32m1_t result1 = __riscv_vfwmul_vf_f32m1( A0, B1, gvl);
+#endif
 
             for(BLASLONG k=1; k < K; ++k) {
                 B0 = B[bi+0];
@@ -565,8 +836,13 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 A0 = __riscv_vle16_v_f16mf2(&A[ai + 0 * gvl], gvl);
                 ai += 4;
 
+#ifdef FP16_NARROW
+                result0 = __riscv_vfmacc_vf_f16mf2(result0, B0, A0, gvl);
+                result1 = __riscv_vfmacc_vf_f16mf2(result1, B1, A0, gvl);
+#else
                 result0 = __riscv_vfwmacc_vf_f32m1(result0, B0, A0, gvl);
                 result1 = __riscv_vfwmacc_vf_f32m1(result1, B1, A0, gvl);
+#endif
             }
 
             BLASLONG ci = n_top * ldc + m_top;
@@ -574,8 +850,13 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat32m1_t c0 = __riscv_vle32_v_f32m1(&C[ci], gvl);
             ci += ldc - gvl * 0;
             vfloat32m1_t c1 = __riscv_vle32_v_f32m1(&C[ci], gvl);
+#ifdef FP16_NARROW
+            c0 = __riscv_vfwmacc_vf_f32m1(c0, alpha, result0, gvl);
+            c1 = __riscv_vfwmacc_vf_f32m1(c1, alpha, result1, gvl);
+#else
             c0 = __riscv_vfmacc_vf_f32m1(c0, alpha, result0, gvl);
             c1 = __riscv_vfmacc_vf_f32m1(c1, alpha, result1, gvl);
+#endif
 
             ci= n_top * ldc + m_top;
 
@@ -591,25 +872,46 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             BLASLONG ai = m_top * K;
             BLASLONG bi = n_top * K;
 
+#ifdef FP16_NARROW
+            _Float16 result0 = 0;
+            _Float16 result1 = 0;
+            _Float16 result2 = 0;
+            _Float16 result3 = 0;
+#else
             float result0 = 0;
             float result1 = 0;
             float result2 = 0;
             float result3 = 0;
+#endif
 
             for(BLASLONG k=0; k<K; k++) {
+#ifdef FP16_NARROW
+                result0+=A[ai+0]*B[bi+0];
+                result1+=A[ai+1]*B[bi+0];
+                result2+=A[ai+0]*B[bi+1];
+                result3+=A[ai+1]*B[bi+1];
+#else
                 result0+=(float)(A[ai+0])*(float)(B[bi+0]);
                 result1+=(float)(A[ai+1])*(float)(B[bi+0]);
                 result2+=(float)(A[ai+0])*(float)(B[bi+1]);
                 result3+=(float)(A[ai+1])*(float)(B[bi+1]);
+#endif
                 ai+=2;
                 bi+=2;
             }
             
             BLASLONG ci=n_top*ldc+m_top;
+#ifdef FP16_NARROW
+            C[ci + 0 * ldc + 0] += alpha * (float)(result0);
+            C[ci + 0 * ldc + 1] += alpha * (float)(result1);
+            C[ci + 1 * ldc + 0] += alpha * (float)(result2);
+            C[ci + 1 * ldc + 1] += alpha * (float)(result3);
+#else
             C[ci + 0 * ldc + 0] += alpha * result0;
             C[ci + 0 * ldc + 1] += alpha * result1;
             C[ci + 1 * ldc + 0] += alpha * result2;
             C[ci + 1 * ldc + 1] += alpha * result3;
+#endif
 
             m_top += 2;
         }
@@ -617,22 +919,37 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
 
         if( M & 1 ) {
             
+#ifdef FP16_NARROW
+            _Float16 result0 = 0;
+            _Float16 result1 = 0;
+#else
             float result0 = 0;
             float result1 = 0;
+#endif
             
             BLASLONG ai = m_top * K;
             BLASLONG bi = n_top * K;
 
             for(BLASLONG k=0; k<K; k++) {
+#ifdef FP16_NARROW
+                result0+=A[ai+0]*B[bi+0];
+                result1+=A[ai+0]*B[bi+1];
+#else
                 result0+=(float)(A[ai+0])*(float)(B[bi+0]);
                 result1+=(float)(A[ai+0])*(float)(B[bi+1]);
+#endif
                 ai+=1;
                 bi+=2;
             }
 
             BLASLONG ci = n_top * ldc + m_top;
+#ifdef FP16_NARROW
+            C[ci + 0 * ldc + 0] += alpha * (float)(result0);
+            C[ci + 1 * ldc + 0] += alpha * (float)(result1);
+#else
             C[ci + 0 * ldc + 0] += alpha * result0;
             C[ci + 1 * ldc + 0] += alpha * result1;
+#endif
             m_top += 1;
         }
 
@@ -656,7 +973,11 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat16m1_t A0 = __riscv_vle16_v_f16m1( &A[ai+0*gvl], gvl );
             ai += 8;
 
+#ifdef FP16_NARROW
+            vfloat16m1_t result0 = __riscv_vfmul_vf_f16m1( A0, B0, gvl);
+#else
             vfloat32m2_t result0 = __riscv_vfwmul_vf_f32m2( A0, B0, gvl);
+#endif
 		
             for(BLASLONG k=1; k<K; k++) {
                 B0 = B[bi+0];
@@ -665,7 +986,11 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 A0 = __riscv_vle16_v_f16m1( &A[ai+0*gvl], gvl );
                 ai += 8;
 
+#ifdef FP16_NARROW
+                result0 = __riscv_vfmacc_vf_f16m1(result0, B0, A0, gvl);
+#else
                 result0 = __riscv_vfwmacc_vf_f32m2(result0, B0, A0, gvl);
+#endif
             }
 
 
@@ -673,7 +998,11 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
 
             vfloat32m2_t c0 = __riscv_vle32_v_f32m2( &C[ci], gvl);
             
+#ifdef FP16_NARROW
+            c0 = __riscv_vfwmacc_vf_f32m2(c0, alpha, result0, gvl);
+#else
             c0 = __riscv_vfmacc_vf_f32m2(c0, alpha, result0, gvl);
+#endif
 
             ci = n_top * ldc + m_top;
 
@@ -692,7 +1021,11 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             vfloat16mf2_t A0 = __riscv_vle16_v_f16mf2(&A[ai + 0 * gvl], gvl);
             ai += 4;
 
+#ifdef FP16_NARROW
+            vfloat16mf2_t result0 = __riscv_vfmul_vf_f16mf2( A0, B0, gvl);
+#else
             vfloat32m1_t result0 = __riscv_vfwmul_vf_f32m1( A0, B0, gvl);
+#endif
 
             for(BLASLONG k=1; k < K; ++k) {
                 B0 = B[bi+0];
@@ -701,13 +1034,21 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
                 A0 = __riscv_vle16_v_f16mf2(&A[ai + 0 * gvl], gvl);
                 ai += 4;
 
+#ifdef FP16_NARROW
+                result0 = __riscv_vfmacc_vf_f16mf2(result0, B0, A0, gvl);
+#else
                 result0 = __riscv_vfwmacc_vf_f32m1(result0, B0, A0, gvl);
+#endif
             }
 
             BLASLONG ci = n_top * ldc + m_top;
 
             vfloat32m1_t c0 = __riscv_vle32_v_f32m1(&C[ci], gvl);
+#ifdef FP16_NARROW
+            c0 = __riscv_vfwmacc_vf_f32m1(c0, alpha, result0, gvl);
+#else
             c0 = __riscv_vfmacc_vf_f32m1(c0, alpha, result0, gvl);
+#endif
 
             ci= n_top * ldc + m_top;
 
@@ -721,20 +1062,35 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
             BLASLONG ai = m_top * K;
             BLASLONG bi = n_top * K;
 
+#ifdef FP16_NARROW
+            _Float16 result0 = 0;
+            _Float16 result1 = 0;
+#else
             float result0 = 0;
             float result1 = 0;
+#endif
 
             for(BLASLONG k=0; k<K; k++) {
+#ifdef FP16_NARROW
+                result0+=A[ai+0]*B[bi+0];
+                result1+=A[ai+1]*B[bi+0];
+#else
                 result0+=(float)(A[ai+0])*(float)(B[bi+0]);
                 result1+=(float)(A[ai+1])*(float)(B[bi+0]);
+#endif
                 ai+=2;
                 bi+=1;
             }
             
             
             BLASLONG ci=n_top*ldc+m_top;
+#ifdef FP16_NARROW
+            C[ci + 0 * ldc + 0] += alpha * (float)(result0);
+            C[ci + 0 * ldc + 1] += alpha * (float)(result1);
+#else
             C[ci + 0 * ldc + 0] += alpha * result0;
             C[ci + 0 * ldc + 1] += alpha * result1;
+#endif
 
             m_top += 2;
         }
@@ -742,19 +1098,31 @@ int CNAME(BLASLONG M, BLASLONG N, BLASLONG K, FLOAT alpha, IFLOAT *A, IFLOAT *B,
 
         if( M & 1 ) {
             
+#ifdef FP16_NARROW
+            _Float16 result0 = 0;
+#else
             float result0 = 0;
+#endif
             
             BLASLONG ai = m_top * K;
             BLASLONG bi = n_top * K;
 
             for(BLASLONG k=0; k<K; k++) {
+#ifdef FP16_NARROW
+                result0+=A[ai+0]*B[bi+0];
+#else
                 result0+=(float)(A[ai+0])*(float)(B[bi+0]);
+#endif
                 ai+=1;
                 bi+=1;
             }
 
             BLASLONG ci = n_top * ldc + m_top;
+#ifdef FP16_NARROW
+            C[ci + 0 * ldc + 0] += alpha * (float)(result0);
+#else
             C[ci + 0 * ldc + 0] += alpha * result0;
+#endif
             m_top += 1;
         }