OSchip
/
OpenBLAS
mirror of https://github.com/xianyi/OpenBLAS.git
0
0
Fork 2
Code Issues Packages Projects Releases Wiki Activity
OpenBLAS/kernel/generic/gemmkernel_4x4.c

265 lines
11 KiB
C
Raw Blame History

/***************************************************************************
* Copyright (c) 2026, 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.
* *****************************************************************************/
/*
* Portable C GEMM micro-kernel with a 4x4 register tile (16 accumulators).
*
* This is a wider companion to gemmkernel_2x2.c intended for in-order scalar
* cores whose FP FMA has a multi-cycle latency but 1/cycle throughput (e.g.
* SiFive U74: fmadd.d latency 7, repeat rate 1). A 2x2 tile exposes only 4
* independent accumulator chains, which is fewer than the FMA latency and
* leaves the FP pipe stalled on the accumulator dependency. A 4x4 tile keeps
* 16 independent chains -- comfortably above the latency -- and lowers the
* load:FMA ratio from 1:1 to 1:2, so the single load/store pipe stops being
* the bottleneck. RV64G has 32 FP registers, so 16 accumulators + 4 A + 4 B
* fit without spilling.
*
* Packed-data contract (identical to the 2x2 kernel, verified against the
* generic tcopy_4 / ncopy_4 copy routines): the A operand is packed by
* tcopy_<UNROLL_M> into MR-row micro-panels [A(r0,k)..A(r3,k)] per k, and the
* B operand by ncopy_<UNROLL_N> into NR-col micro-panels [B(k,c0)..B(k,c3)]
* per k. Both dimensions are decomposed as 4 / 2 / 1 sub-blocks at the edges.
*/
#include "common.h"
#include "conversion_macros.h"
#ifdef BGEMM
#define C_TO_F32 TO_F32
#else
#define C_TO_F32
#endif
int CNAME(BLASLONG bm,BLASLONG bn,BLASLONG bk,FLOAT alpha,IFLOAT* ba,IFLOAT* bb,FLOAT* C,BLASLONG ldc
#ifdef TRMMKERNEL
,BLASLONG offset
#endif
)
{
BLASLONG i,j,k;
FLOAT *C0,*C1,*C2,*C3;
IFLOAT *ptrba,*ptrbb;
FLOAT r0c0,r1c0,r2c0,r3c0;
FLOAT r0c1,r1c1,r2c1,r3c1;
FLOAT r0c2,r1c2,r2c2,r3c2;
FLOAT r0c3,r1c3,r2c3,r3c3;
IFLOAT a0,a1,a2,a3,b0,b1,b2,b3;
/* ==================== N panels of 4 ==================== */
for (j=0; j<bn/4; j+=1)
{
C0 = C;
C1 = C0+ldc;
C2 = C1+ldc;
C3 = C2+ldc;
ptrba = ba;
/* ---- 4x4 : 4 rows x 4 cols, 16 accumulators ---- */
for (i=0; i<bm/4; i+=1)
{
ptrbb = bb;
r0c0=r1c0=r2c0=r3c0=0;
r0c1=r1c1=r2c1=r3c1=0;
r0c2=r1c2=r2c2=r3c2=0;
r0c3=r1c3=r2c3=r3c3=0;
for (k=0; k<bk; k+=1)
{
b0=ptrbb[0]; b1=ptrbb[1]; b2=ptrbb[2]; b3=ptrbb[3];
a0=ptrba[0]; a1=ptrba[1]; a2=ptrba[2]; a3=ptrba[3];
r0c0+=TO_F32(a0)*TO_F32(b0); r1c0+=TO_F32(a1)*TO_F32(b0); r2c0+=TO_F32(a2)*TO_F32(b0); r3c0+=TO_F32(a3)*TO_F32(b0);
r0c1+=TO_F32(a0)*TO_F32(b1); r1c1+=TO_F32(a1)*TO_F32(b1); r2c1+=TO_F32(a2)*TO_F32(b1); r3c1+=TO_F32(a3)*TO_F32(b1);
r0c2+=TO_F32(a0)*TO_F32(b2); r1c2+=TO_F32(a1)*TO_F32(b2); r2c2+=TO_F32(a2)*TO_F32(b2); r3c2+=TO_F32(a3)*TO_F32(b2);
r0c3+=TO_F32(a0)*TO_F32(b3); r1c3+=TO_F32(a1)*TO_F32(b3); r2c3+=TO_F32(a2)*TO_F32(b3); r3c3+=TO_F32(a3)*TO_F32(b3);
ptrba+=4; ptrbb+=4;
}
C0[0]=TO_OUTPUT(C_TO_F32(C0[0])+r0c0*ALPHA); C0[1]=TO_OUTPUT(C_TO_F32(C0[1])+r1c0*ALPHA); C0[2]=TO_OUTPUT(C_TO_F32(C0[2])+r2c0*ALPHA); C0[3]=TO_OUTPUT(C_TO_F32(C0[3])+r3c0*ALPHA);
C1[0]=TO_OUTPUT(C_TO_F32(C1[0])+r0c1*ALPHA); C1[1]=TO_OUTPUT(C_TO_F32(C1[1])+r1c1*ALPHA); C1[2]=TO_OUTPUT(C_TO_F32(C1[2])+r2c1*ALPHA); C1[3]=TO_OUTPUT(C_TO_F32(C1[3])+r3c1*ALPHA);
C2[0]=TO_OUTPUT(C_TO_F32(C2[0])+r0c2*ALPHA); C2[1]=TO_OUTPUT(C_TO_F32(C2[1])+r1c2*ALPHA); C2[2]=TO_OUTPUT(C_TO_F32(C2[2])+r2c2*ALPHA); C2[3]=TO_OUTPUT(C_TO_F32(C2[3])+r3c2*ALPHA);
C3[0]=TO_OUTPUT(C_TO_F32(C3[0])+r0c3*ALPHA); C3[1]=TO_OUTPUT(C_TO_F32(C3[1])+r1c3*ALPHA); C3[2]=TO_OUTPUT(C_TO_F32(C3[2])+r2c3*ALPHA); C3[3]=TO_OUTPUT(C_TO_F32(C3[3])+r3c3*ALPHA);
C0+=4; C1+=4; C2+=4; C3+=4;
}
/* ---- 2x4 : 2 rows x 4 cols ---- */
if (bm & 2)
{
ptrbb = bb;
r0c0=r1c0=0; r0c1=r1c1=0; r0c2=r1c2=0; r0c3=r1c3=0;
for (k=0; k<bk; k+=1)
{
b0=ptrbb[0]; b1=ptrbb[1]; b2=ptrbb[2]; b3=ptrbb[3];
a0=ptrba[0]; a1=ptrba[1];
r0c0+=TO_F32(a0)*TO_F32(b0); r1c0+=TO_F32(a1)*TO_F32(b0);
r0c1+=TO_F32(a0)*TO_F32(b1); r1c1+=TO_F32(a1)*TO_F32(b1);
r0c2+=TO_F32(a0)*TO_F32(b2); r1c2+=TO_F32(a1)*TO_F32(b2);
r0c3+=TO_F32(a0)*TO_F32(b3); r1c3+=TO_F32(a1)*TO_F32(b3);
ptrba+=2; ptrbb+=4;
}
C0[0]=TO_OUTPUT(C_TO_F32(C0[0])+r0c0*ALPHA); C0[1]=TO_OUTPUT(C_TO_F32(C0[1])+r1c0*ALPHA);
C1[0]=TO_OUTPUT(C_TO_F32(C1[0])+r0c1*ALPHA); C1[1]=TO_OUTPUT(C_TO_F32(C1[1])+r1c1*ALPHA);
C2[0]=TO_OUTPUT(C_TO_F32(C2[0])+r0c2*ALPHA); C2[1]=TO_OUTPUT(C_TO_F32(C2[1])+r1c2*ALPHA);
C3[0]=TO_OUTPUT(C_TO_F32(C3[0])+r0c3*ALPHA); C3[1]=TO_OUTPUT(C_TO_F32(C3[1])+r1c3*ALPHA);
C0+=2; C1+=2; C2+=2; C3+=2;
}
/* ---- 1x4 : 1 row x 4 cols ---- */
if (bm & 1)
{
ptrbb = bb;
r0c0=0; r0c1=0; r0c2=0; r0c3=0;
for (k=0; k<bk; k+=1)
{
b0=ptrbb[0]; b1=ptrbb[1]; b2=ptrbb[2]; b3=ptrbb[3];
a0=ptrba[0];
r0c0+=TO_F32(a0)*TO_F32(b0);
r0c1+=TO_F32(a0)*TO_F32(b1);
r0c2+=TO_F32(a0)*TO_F32(b2);
r0c3+=TO_F32(a0)*TO_F32(b3);
ptrba+=1; ptrbb+=4;
}
C0[0]=TO_OUTPUT(C_TO_F32(C0[0])+r0c0*ALPHA);
C1[0]=TO_OUTPUT(C_TO_F32(C1[0])+r0c1*ALPHA);
C2[0]=TO_OUTPUT(C_TO_F32(C2[0])+r0c2*ALPHA);
C3[0]=TO_OUTPUT(C_TO_F32(C3[0])+r0c3*ALPHA);
C0+=1; C1+=1; C2+=1; C3+=1;
}
bb = bb + bk*4;
C = C + ldc*4;
}
/* ==================== N panel of 2 ==================== */
if (bn & 2)
{
C0 = C;
C1 = C0+ldc;
ptrba = ba;
for (i=0; i<bm/4; i+=1)
{
ptrbb = bb;
r0c0=r1c0=r2c0=r3c0=0;
r0c1=r1c1=r2c1=r3c1=0;
for (k=0; k<bk; k+=1)
{
b0=ptrbb[0]; b1=ptrbb[1];
a0=ptrba[0]; a1=ptrba[1]; a2=ptrba[2]; a3=ptrba[3];
r0c0+=TO_F32(a0)*TO_F32(b0); r1c0+=TO_F32(a1)*TO_F32(b0); r2c0+=TO_F32(a2)*TO_F32(b0); r3c0+=TO_F32(a3)*TO_F32(b0);
r0c1+=TO_F32(a0)*TO_F32(b1); r1c1+=TO_F32(a1)*TO_F32(b1); r2c1+=TO_F32(a2)*TO_F32(b1); r3c1+=TO_F32(a3)*TO_F32(b1);
ptrba+=4; ptrbb+=2;
}
C0[0]=TO_OUTPUT(C_TO_F32(C0[0])+r0c0*ALPHA); C0[1]=TO_OUTPUT(C_TO_F32(C0[1])+r1c0*ALPHA); C0[2]=TO_OUTPUT(C_TO_F32(C0[2])+r2c0*ALPHA); C0[3]=TO_OUTPUT(C_TO_F32(C0[3])+r3c0*ALPHA);
C1[0]=TO_OUTPUT(C_TO_F32(C1[0])+r0c1*ALPHA); C1[1]=TO_OUTPUT(C_TO_F32(C1[1])+r1c1*ALPHA); C1[2]=TO_OUTPUT(C_TO_F32(C1[2])+r2c1*ALPHA); C1[3]=TO_OUTPUT(C_TO_F32(C1[3])+r3c1*ALPHA);
C0+=4; C1+=4;
}
if (bm & 2)
{
ptrbb = bb;
r0c0=r1c0=0; r0c1=r1c1=0;
for (k=0; k<bk; k+=1)
{
b0=ptrbb[0]; b1=ptrbb[1];
a0=ptrba[0]; a1=ptrba[1];
r0c0+=TO_F32(a0)*TO_F32(b0); r1c0+=TO_F32(a1)*TO_F32(b0);
r0c1+=TO_F32(a0)*TO_F32(b1); r1c1+=TO_F32(a1)*TO_F32(b1);
ptrba+=2; ptrbb+=2;
}
C0[0]=TO_OUTPUT(C_TO_F32(C0[0])+r0c0*ALPHA); C0[1]=TO_OUTPUT(C_TO_F32(C0[1])+r1c0*ALPHA);
C1[0]=TO_OUTPUT(C_TO_F32(C1[0])+r0c1*ALPHA); C1[1]=TO_OUTPUT(C_TO_F32(C1[1])+r1c1*ALPHA);
C0+=2; C1+=2;
}
if (bm & 1)
{
ptrbb = bb;
r0c0=0; r0c1=0;
for (k=0; k<bk; k+=1)
{
b0=ptrbb[0]; b1=ptrbb[1];
a0=ptrba[0];
r0c0+=TO_F32(a0)*TO_F32(b0);
r0c1+=TO_F32(a0)*TO_F32(b1);
ptrba+=1; ptrbb+=2;
}
C0[0]=TO_OUTPUT(C_TO_F32(C0[0])+r0c0*ALPHA);
C1[0]=TO_OUTPUT(C_TO_F32(C1[0])+r0c1*ALPHA);
C0+=1; C1+=1;
}
bb = bb + bk*2;
C = C + ldc*2;
}
/* ==================== N panel of 1 ==================== */
if (bn & 1)
{
C0 = C;
ptrba = ba;
for (i=0; i<bm/4; i+=1)
{
ptrbb = bb;
r0c0=r1c0=r2c0=r3c0=0;
for (k=0; k<bk; k+=1)
{
b0=ptrbb[0];
a0=ptrba[0]; a1=ptrba[1]; a2=ptrba[2]; a3=ptrba[3];
r0c0+=TO_F32(a0)*TO_F32(b0); r1c0+=TO_F32(a1)*TO_F32(b0); r2c0+=TO_F32(a2)*TO_F32(b0); r3c0+=TO_F32(a3)*TO_F32(b0);
ptrba+=4; ptrbb+=1;
}
C0[0]=TO_OUTPUT(C_TO_F32(C0[0])+r0c0*ALPHA); C0[1]=TO_OUTPUT(C_TO_F32(C0[1])+r1c0*ALPHA); C0[2]=TO_OUTPUT(C_TO_F32(C0[2])+r2c0*ALPHA); C0[3]=TO_OUTPUT(C_TO_F32(C0[3])+r3c0*ALPHA);
C0+=4;
}
if (bm & 2)
{
ptrbb = bb;
r0c0=r1c0=0;
for (k=0; k<bk; k+=1)
{
b0=ptrbb[0];
a0=ptrba[0]; a1=ptrba[1];
r0c0+=TO_F32(a0)*TO_F32(b0); r1c0+=TO_F32(a1)*TO_F32(b0);
ptrba+=2; ptrbb+=1;
}
C0[0]=TO_OUTPUT(C_TO_F32(C0[0])+r0c0*ALPHA); C0[1]=TO_OUTPUT(C_TO_F32(C0[1])+r1c0*ALPHA);
C0+=2;
}
if (bm & 1)
{
ptrbb = bb;
r0c0=0;
for (k=0; k<bk; k+=1)
{
r0c0+=TO_F32(ptrba[0])*TO_F32(ptrbb[0]);
ptrba+=1; ptrbb+=1;
}
C0[0]=TO_OUTPUT(C_TO_F32(C0[0])+r0c0*ALPHA);
C0+=1;
}
bb = bb + bk;
C = C + ldc;
}
return 0;
}
Reference in New Issue View Git Blame Copy Permalink