youmuliu/stgsharp
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Double add sub mul sc-mul ok

Browse Source
This commit is contained in:
Youmu Liu
2025-12-31 18:15:43 +08:00
parent 6484307b8a
commit 0e93787fa2
18 changed files with 1210 additions and 944 deletions
Download Patch File Download Diff File Expand all files Collapse all files
StgSharp.Native/sn-head/sn_intrinsic.h
+5
View File
@@ -264,6 +264,11 @@ DECLARE_PNL_PROC_LEFT_RIGHT_ANS(float, avx, fma);
DECLARE_PNL_PROC_LEFT_RIGHT_ANS(float, avx_fma, fma);
DECLARE_PNL_PROC_LEFT_RIGHT_ANS(float, 512, fma);
DECLARE_PNL_PROC_LEFT_RIGHT_ANS(double, sse, fma);
DECLARE_PNL_PROC_LEFT_RIGHT_ANS(double, avx, fma);
DECLARE_PNL_PROC_LEFT_RIGHT_ANS(double, avx_fma, fma);
DECLARE_PNL_PROC_LEFT_RIGHT_ANS(double, 512, fma);
#pragma endregion
#pragma region scaler_simd
StgSharp.Native/sn-src/k_fma.c
+177 -1
View File
@@ -142,7 +142,7 @@ DECLARE_PNL_PROC_LEFT_RIGHT_ANS(float, 512, fma)
{
const __m512i ones = _mm512_set1_epi32(1);
for (int col = 0; col < 8; ++col) {
for (int col = 0; col < 16; ++col) {
register __m512 bcol = _mm512_load_ps((float *)&right->vec[col][0]);
register __m512 ccol = _mm512_load_ps((float *)&ans->vec[col][0]);
register __m512i idx = _mm512_setzero_si512();
@@ -203,3 +203,179 @@ DECLARE_PNL_PROC_LEFT_RIGHT_ANS(float, 512, fma)
_mm512_store_ps((float *)&ans->vec[col][0], ccol);
}
}
DECLARE_PNL_PROC_LEFT_RIGHT_ANS(double, sse, fma)
{
for (int col = 0; col < 4; ++col) {
/* current answer column (two chunks: rows 0-1 and 2-3) */
register __m128d c_lo = _mm_castps_pd(PANEL_COL_VEC((*ans), col, 0));
register __m128d c_hi = _mm_castps_pd(PANEL_COL_VEC((*ans), col, 1));
/* right matrix column split into two chunks */
register __m128d b_lo = _mm_castps_pd(PANEL_COL_VEC((*right), col, 0));
register __m128d b_hi = _mm_castps_pd(PANEL_COL_VEC((*right), col, 1));
/* k = 0 */
{
register __m128d a_lo = _mm_castps_pd(PANEL_COL_VEC((*left), 0, 0));
register __m128d a_hi = _mm_castps_pd(PANEL_COL_VEC((*left), 0, 1));
register __m128d b = _mm_shuffle_pd(b_lo, b_lo, 0x0);
c_lo = _mm_add_pd(c_lo, _mm_mul_pd(a_lo, b));
c_hi = _mm_add_pd(c_hi, _mm_mul_pd(a_hi, b));
}
/* k = 1 */
{
register __m128d a_lo = _mm_castps_pd(PANEL_COL_VEC((*left), 1, 0));
register __m128d a_hi = _mm_castps_pd(PANEL_COL_VEC((*left), 1, 1));
register __m128d b = _mm_shuffle_pd(b_lo, b_lo, 0x3);
c_lo = _mm_add_pd(c_lo, _mm_mul_pd(a_lo, b));
c_hi = _mm_add_pd(c_hi, _mm_mul_pd(a_hi, b));
}
/* k = 2 */
{
register __m128d a_lo = _mm_castps_pd(PANEL_COL_VEC((*left), 2, 0));
register __m128d a_hi = _mm_castps_pd(PANEL_COL_VEC((*left), 2, 1));
register __m128d b = _mm_shuffle_pd(b_hi, b_hi, 0x0);
c_lo = _mm_add_pd(c_lo, _mm_mul_pd(a_lo, b));
c_hi = _mm_add_pd(c_hi, _mm_mul_pd(a_hi, b));
}
/* k = 3 */
{
register __m128d a_lo = _mm_castps_pd(PANEL_COL_VEC((*left), 3, 0));
register __m128d a_hi = _mm_castps_pd(PANEL_COL_VEC((*left), 3, 1));
register __m128d b = _mm_shuffle_pd(b_hi, b_hi, 0x3);
c_lo = _mm_add_pd(c_lo, _mm_mul_pd(a_lo, b));
c_hi = _mm_add_pd(c_hi, _mm_mul_pd(a_hi, b));
}
PANEL_COL_VEC((*ans), col, 0) = _mm_castpd_ps(c_lo);
PANEL_COL_VEC((*ans), col, 1) = _mm_castpd_ps(c_hi);
}
}
DECLARE_PNL_PROC_LEFT_RIGHT_ANS(double, avx, fma)
{
register __m256d a0 = _mm256_load_pd((double *)&left->vec[0][0]);
register __m256d a1 = _mm256_load_pd((double *)&left->vec[1][0]);
register __m256d a2 = _mm256_load_pd((double *)&left->vec[2][0]);
register __m256d a3 = _mm256_load_pd((double *)&left->vec[3][0]);
register __m256d c0 = _mm256_load_pd((double *)&ans->vec[0][0]);
register __m256d c1 = _mm256_load_pd((double *)&ans->vec[1][0]);
register __m256d c2 = _mm256_load_pd((double *)&ans->vec[2][0]);
register __m256d c3 = _mm256_load_pd((double *)&ans->vec[3][0]);
register __m256d bcol = _mm256_load_pd((double *)&right->vec[0][0]);
c0 = _mm256_add_pd(c0, _mm256_mul_pd(a0, _mm256_permute_pd(bcol, 0x0)));
c0 = _mm256_add_pd(c0, _mm256_mul_pd(a1, _mm256_permute_pd(bcol, 0x3)));
c0 = _mm256_add_pd(c0, _mm256_mul_pd(a2, _mm256_permute_pd(bcol, 0x0)));
c0 = _mm256_add_pd(c0, _mm256_mul_pd(a3, _mm256_permute_pd(bcol, 0x3)));
_mm256_store_pd((double *)&ans->vec[0][0], c0);
bcol = _mm256_load_pd((double *)&right->vec[0][0]);
c1 = _mm256_add_pd(c1, _mm256_mul_pd(a0, _mm256_permute_pd(bcol, 0x0)));
c1 = _mm256_add_pd(c1, _mm256_mul_pd(a1, _mm256_permute_pd(bcol, 0x3)));
c1 = _mm256_add_pd(c1, _mm256_mul_pd(a2, _mm256_permute_pd(bcol, 0x0)));
c1 = _mm256_add_pd(c1, _mm256_mul_pd(a3, _mm256_permute_pd(bcol, 0x3)));
_mm256_store_pd((double *)&ans->vec[1][0], c1);
bcol = _mm256_load_pd((double *)&right->vec[0][0]);
c2 = _mm256_add_pd(c2, _mm256_mul_pd(a0, _mm256_permute_pd(bcol, 0x0)));
c2 = _mm256_add_pd(c2, _mm256_mul_pd(a1, _mm256_permute_pd(bcol, 0x3)));
c2 = _mm256_add_pd(c2, _mm256_mul_pd(a2, _mm256_permute_pd(bcol, 0x0)));
c2 = _mm256_add_pd(c2, _mm256_mul_pd(a3, _mm256_permute_pd(bcol, 0x3)));
_mm256_store_pd((double *)&ans->vec[2][0], c2);
bcol = _mm256_load_pd((double *)&right->vec[0][0]);
c3 = _mm256_add_pd(c3, _mm256_mul_pd(a0, _mm256_permute_pd(bcol, 0x0)));
c3 = _mm256_add_pd(c3, _mm256_mul_pd(a1, _mm256_permute_pd(bcol, 0x3)));
c3 = _mm256_add_pd(c3, _mm256_mul_pd(a2, _mm256_permute_pd(bcol, 0x0)));
c3 = _mm256_add_pd(c3, _mm256_mul_pd(a3, _mm256_permute_pd(bcol, 0x3)));
_mm256_store_pd((double *)&ans->vec[3][0], c3);
}
DECLARE_PNL_PROC_LEFT_RIGHT_ANS(double, avx_fma, fma)
{
register __m256d a0 = _mm256_load_pd((double *)&left->vec[0][0]);
register __m256d a1 = _mm256_load_pd((double *)&left->vec[1][0]);
register __m256d a2 = _mm256_load_pd((double *)&left->vec[2][0]);
register __m256d a3 = _mm256_load_pd((double *)&left->vec[3][0]);
register __m256d c0 = _mm256_load_pd((double *)&ans->vec[0][0]);
register __m256d c1 = _mm256_load_pd((double *)&ans->vec[1][0]);
register __m256d c2 = _mm256_load_pd((double *)&ans->vec[2][0]);
register __m256d c3 = _mm256_load_pd((double *)&ans->vec[3][0]);
register __m256d bcol = _mm256_load_pd((double *)&right->vec[0][0]);
c0 = _mm256_fmadd_pd(a0, _mm256_permute_pd(bcol, 0x0), c0);
c0 = _mm256_fmadd_pd(a1, _mm256_permute_pd(bcol, 0x3), c0);
c0 = _mm256_fmadd_pd(a2, _mm256_permute_pd(bcol, 0x0), c0);
c0 = _mm256_fmadd_pd(a3, _mm256_permute_pd(bcol, 0x3), c0);
_mm256_store_pd((double *)&ans->vec[0][0], c0);
bcol = _mm256_load_pd((double *)&right->vec[0][0]);
c1 = _mm256_fmadd_pd(a0, _mm256_permute_pd(bcol, 0x0), c1);
c1 = _mm256_fmadd_pd(a1, _mm256_permute_pd(bcol, 0x3), c1);
c1 = _mm256_fmadd_pd(a2, _mm256_permute_pd(bcol, 0x0), c1);
c1 = _mm256_fmadd_pd(a3, _mm256_permute_pd(bcol, 0x3), c1);
_mm256_store_pd((double *)&ans->vec[1][0], c1);
bcol = _mm256_load_pd((double *)&right->vec[0][0]);
c2 = _mm256_fmadd_pd(a0, _mm256_permute_pd(bcol, 0x0), c2);
c2 = _mm256_fmadd_pd(a1, _mm256_permute_pd(bcol, 0x3), c2);
c2 = _mm256_fmadd_pd(a2, _mm256_permute_pd(bcol, 0x0), c2);
c2 = _mm256_fmadd_pd(a3, _mm256_permute_pd(bcol, 0x3), c2);
_mm256_store_pd((double *)&ans->vec[2][0], c2);
bcol = _mm256_load_pd((double *)&right->vec[0][0]);
c3 = _mm256_fmadd_pd(a0, _mm256_permute_pd(bcol, 0x0), c3);
c3 = _mm256_fmadd_pd(a1, _mm256_permute_pd(bcol, 0x3), c3);
c3 = _mm256_fmadd_pd(a2, _mm256_permute_pd(bcol, 0x0), c3);
c3 = _mm256_fmadd_pd(a3, _mm256_permute_pd(bcol, 0x3), c3);
_mm256_store_pd((double *)&ans->vec[3][0], c3);
}
// ans = left * right + ans
// right is not transposed
DECLARE_PNL_PROC_LEFT_RIGHT_ANS(double, 512, fma)
{
const __m512i ones = _mm512_set1_epi32(1);
register __m512d a0 = _mm512_load_pd((double *)&left->vec[0][0]);
register __m512d a1 = _mm512_load_pd((double *)&left->vec[1][0]);
register __m512d a2 = _mm512_load_pd((double *)&left->vec[2][0]);
register __m512d a3 = _mm512_load_pd((double *)&left->vec[3][0]);
register __m512d a4 = _mm512_load_pd((double *)&left->vec[4][0]);
register __m512d a5 = _mm512_load_pd((double *)&left->vec[5][0]);
register __m512d a6 = _mm512_load_pd((double *)&left->vec[6][0]);
register __m512d a7 = _mm512_load_pd((double *)&left->vec[7][0]);
for (int col = 0; col < 8; ++col) {
register __m512d bcol = _mm512_load_pd((double *)&right->vec[col][0]);
register __m512d ccol = _mm512_load_pd((double *)&ans->vec[col][0]);
register __m512i idx = _mm512_setzero_si512();
register __m512i increament = _mm512_load_si512(&ones);
ccol = _mm512_fmadd_pd(a0, _mm512_permutexvar_pd(idx, bcol), ccol);
idx = _mm512_add_epi64(idx, increament);
ccol = _mm512_fmadd_pd(a1, _mm512_permutexvar_pd(idx, bcol), ccol);
idx = _mm512_add_epi64(idx, increament);
ccol = _mm512_fmadd_pd(a2, _mm512_permutexvar_pd(idx, bcol), ccol);
idx = _mm512_add_epi64(idx, increament);
ccol = _mm512_fmadd_pd(a3, _mm512_permutexvar_pd(idx, bcol), ccol);
idx = _mm512_add_epi64(idx, increament);
ccol = _mm512_fmadd_pd(a4, _mm512_permutexvar_pd(idx, bcol), ccol);
idx = _mm512_add_epi64(idx, increament);
ccol = _mm512_fmadd_pd(a5, _mm512_permutexvar_pd(idx, bcol), ccol);
idx = _mm512_add_epi64(idx, increament);
ccol = _mm512_fmadd_pd(a6, _mm512_permutexvar_pd(idx, bcol), ccol);
idx = _mm512_add_epi64(idx, increament);
ccol = _mm512_fmadd_pd(a7, _mm512_permutexvar_pd(idx, bcol), ccol);
idx = _mm512_add_epi64(idx, increament);
_mm512_store_pd((double *)&ans->vec[col][0], ccol);
}
}
StgSharp.Native/sn-src/load_intrinsic.c
+3 -3
View File
@@ -84,7 +84,7 @@ void load_intrinsic_sse(sn_intrinsic *context)
context->mat[F64].buffer_transpose = NULL;
context->mat[F64].build_panel = NULL;
context->mat[F64].clear_panel = PNL_PROC(double, sse, clear);
context->mat[F64].panel_fma = NULL;
context->mat[F64].panel_fma = PNL_PROC(double, sse, fma);
context->mat[F64].pivot = NULL;
context->mat[F64].store_panel = NULL;
@@ -115,7 +115,7 @@ void load_intrinsic_avx(sn_intrinsic *context)
context->mat[F64].buffer_transpose = NULL;
context->mat[F64].build_panel = NULL;
context->mat[F64].clear_panel = PNL_PROC(double, avx, clear);
context->mat[F64].panel_fma = NULL;
context->mat[F64].panel_fma = PNL_PROC(double, avx, fma);
context->mat[F64].pivot = NULL;
context->mat[F64].store_panel = NULL;
@@ -146,7 +146,7 @@ void load_intrinsic_512(sn_intrinsic *context)
context->mat[F64].buffer_transpose = NULL;
context->mat[F64].build_panel = BUILD_PANEL(512, double);
context->mat[F64].clear_panel = PNL_PROC(double, 512, clear);
context->mat[F64].panel_fma = NULL;
context->mat[F64].panel_fma = PNL_PROC(double, 512, fma);
context->mat[F64].pivot = NULL;
context->mat[F64].store_panel = STORE_PANEL(512, double);