atl_smmmncu_av.c

来自「基于Blas CLapck的.用过的人知道是干啥的」· C语言 代码 · 共 484 行 · 第 1/2 页

         vB3 = vec_ld(0, pB0+KB3); pB0 += 4;

         for (k=kstart; k; k--)
         {
            vC0_0 = vec_madd(vA0, vB0, vC0_0);
            vC1_0 = vec_madd(vA1, vB0, vC1_0);
            vC2_0 = vec_madd(vA2, vB0, vC2_0);
            vC3_0 = vec_madd(vA3, vB0, vC3_0); vB0 = vec_ld(0, pB0);
            vC0_1 = vec_madd(vA0, vB1, vC0_1);
            vC1_1 = vec_madd(vA1, vB1, vC1_1);
            vC2_1 = vec_madd(vA2, vB1, vC2_1);
            vC3_1 = vec_madd(vA3, vB1, vC3_1); vB1 = vec_ld(0, pB0+KB);
            vC0_2 = vec_madd(vA0, vB2, vC0_2);
            vC1_2 = vec_madd(vA1, vB2, vC1_2);
            vC2_2 = vec_madd(vA2, vB2, vC2_2);
            vC3_2 = vec_madd(vA3, vB2, vC3_2); vB2 = vec_ld(0, pB0+KB2);
            vC0_3 = vec_madd(vA0, vB3, vC0_3); vA0 = vec_ld(0, pA0);
            vC1_3 = vec_madd(vA1, vB3, vC1_3); vA1 = vec_ld(0, pA0+KB);
            vC2_3 = vec_madd(vA2, vB3, vC2_3); vA2 = vec_ld(0, pA0+KB2);
            vC3_3 = vec_madd(vA3, vB3, vC3_3); vA3 = vec_ld(0, pA0+KB3);

            vC0_0 = vec_madd(vA0, vB0, vC0_0); vB3 = vec_ld(0, pB0+KB3);
            vC1_0 = vec_madd(vA1, vB0, vC1_0);
            vC2_0 = vec_madd(vA2, vB0, vC2_0);
            vC3_0 = vec_madd(vA3, vB0, vC3_0); vB0 = vec_ld(0, pB0+4);
            vC0_1 = vec_madd(vA0, vB1, vC0_1);
            vC1_1 = vec_madd(vA1, vB1, vC1_1);
            vC2_1 = vec_madd(vA2, vB1, vC2_1);
            vC3_1 = vec_madd(vA3, vB1, vC3_1); vB1 = vec_ld(0, pB0+KB+4);
            vC0_2 = vec_madd(vA0, vB2, vC0_2);
            vC1_2 = vec_madd(vA1, vB2, vC1_2);
            vC2_2 = vec_madd(vA2, vB2, vC2_2);
            vC3_2 = vec_madd(vA3, vB2, vC3_2); vB2 = vec_ld(0, pB0+KB2+4);
            vC0_3 = vec_madd(vA0, vB3, vC0_3); vA0 = vec_ld(0, pA0+4);
            vC1_3 = vec_madd(vA1, vB3, vC1_3); vA1 = vec_ld(0, pA0+KB+4);
            vC2_3 = vec_madd(vA2, vB3, vC2_3); vA2 = vec_ld(0, pA0+KB2+4);
            vC3_3 = vec_madd(vA3, vB3, vC3_3); vA3 = vec_ld(0, pA0+KB3+4);

            vC0_0 = vec_madd(vA0, vB0, vC0_0); vB3 = vec_ld(0, pB0+KB3+4);
            vC1_0 = vec_madd(vA1, vB0, vC1_0);
            vC2_0 = vec_madd(vA2, vB0, vC2_0);
            vC3_0 = vec_madd(vA3, vB0, vC3_0); vB0 = vec_ld(0, pB0+8);
            vC0_1 = vec_madd(vA0, vB1, vC0_1);
            vC1_1 = vec_madd(vA1, vB1, vC1_1);
            vC2_1 = vec_madd(vA2, vB1, vC2_1);
            vC3_1 = vec_madd(vA3, vB1, vC3_1); vB1 = vec_ld(0, pB0+KB+8);
            vC0_2 = vec_madd(vA0, vB2, vC0_2);
            vC1_2 = vec_madd(vA1, vB2, vC1_2);
            vC2_2 = vec_madd(vA2, vB2, vC2_2);
            vC3_2 = vec_madd(vA3, vB2, vC3_2); vB2 = vec_ld(0, pB0+KB2+8);
            vC0_3 = vec_madd(vA0, vB3, vC0_3); vA0 = vec_ld(0, pA0+8);
            vC1_3 = vec_madd(vA1, vB3, vC1_3); vA1 = vec_ld(0, pA0+KB+8);
            vC2_3 = vec_madd(vA2, vB3, vC2_3); vA2 = vec_ld(0, pA0+KB2+8);
            vC3_3 = vec_madd(vA3, vB3, vC3_3); vA3 = vec_ld(0, pA0+KB3+8);

            vC0_0 = vec_madd(vA0, vB0, vC0_0); vB3 = vec_ld(0, pB0+KB3+8);
            vC1_0 = vec_madd(vA1, vB0, vC1_0);
            vC2_0 = vec_madd(vA2, vB0, vC2_0);
            vC3_0 = vec_madd(vA3, vB0, vC3_0); vB0 = vec_ld(0, pB0+12);
            vC0_1 = vec_madd(vA0, vB1, vC0_1);
            vC1_1 = vec_madd(vA1, vB1, vC1_1);
            vC2_1 = vec_madd(vA2, vB1, vC2_1);
            vC3_1 = vec_madd(vA3, vB1, vC3_1); vB1 = vec_ld(0, pB0+KB+12);
            vC0_2 = vec_madd(vA0, vB2, vC0_2);
            vC1_2 = vec_madd(vA1, vB2, vC1_2);
            vC2_2 = vec_madd(vA2, vB2, vC2_2);
            vC3_2 = vec_madd(vA3, vB2, vC3_2); vB2 = vec_ld(0, pB0+KB2+12);
            vC0_3 = vec_madd(vA0, vB3, vC0_3); vA0 = vec_ld(0, pA0+12);
            vC1_3 = vec_madd(vA1, vB3, vC1_3); vA1 = vec_ld(0, pA0+KB+12);
            vC2_3 = vec_madd(vA2, vB3, vC2_3); vA2 = vec_ld(0, pA0+KB2+12);
            vC3_3 = vec_madd(vA3, vB3, vC3_3); vA3 = vec_ld(0, pA0+KB3+12);
                  vB3 = vec_ld(0, pB0+KB3+12);
            pA0 += 16;
            pB0 += 16;
         }
         vC0_0 = vec_madd(vA0, vB0, vC0_0); ATL_pfavW(pC0, cwrdC, 0);
         vC1_0 = vec_madd(vA1, vB0, vC1_0); ATL_pfavW(pC1, cwrdC, 1);
         vC2_0 = vec_madd(vA2, vB0, vC2_0); ATL_pfavW(pC2, cwrdC, 2);
         vC3_0 = vec_madd(vA3, vB0, vC3_0); vB0 = vec_ld(0, pB0);
         vC0_1 = vec_madd(vA0, vB1, vC0_1); ATL_pfavW(pC3, cwrdC, 3);
         vC1_1 = vec_madd(vA1, vB1, vC1_1);
         vC2_1 = vec_madd(vA2, vB1, vC2_1);
         vC3_1 = vec_madd(vA3, vB1, vC3_1); vB1 = vec_ld(0, pB0+KB);
         vC0_2 = vec_madd(vA0, vB2, vC0_2);
         vC1_2 = vec_madd(vA1, vB2, vC1_2);
         vC2_2 = vec_madd(vA2, vB2, vC2_2);
         vC3_2 = vec_madd(vA3, vB2, vC3_2); vB2 = vec_ld(0, pB0+KB2);
         vC0_3 = vec_madd(vA0, vB3, vC0_3); vA0 = vec_ld(0, pA0);
         vC1_3 = vec_madd(vA1, vB3, vC1_3); vA1 = vec_ld(0, pA0+KB);
         vC2_3 = vec_madd(vA2, vB3, vC2_3); vA2 = vec_ld(0, pA0+KB2);
         vC3_3 = vec_madd(vA3, vB3, vC3_3); vA3 = vec_ld(0, pA0+KB3);

         vC0_0 = vec_madd(vA0, vB0, vC0_0); vB3 = vec_ld(0, pB0+KB3);
         vC1_0 = vec_madd(vA1, vB0, vC1_0);
         vC2_0 = vec_madd(vA2, vB0, vC2_0);
         vC3_0 = vec_madd(vA3, vB0, vC3_0); vB0 = vec_ld(0, pB0+4);
         vC0_1 = vec_madd(vA0, vB1, vC0_1);
         vC1_1 = vec_madd(vA1, vB1, vC1_1);
         vC2_1 = vec_madd(vA2, vB1, vC2_1);
         vC3_1 = vec_madd(vA3, vB1, vC3_1); vB1 = vec_ld(0, pB0+KB+4);
         vC0_2 = vec_madd(vA0, vB2, vC0_2);
         vC1_2 = vec_madd(vA1, vB2, vC1_2);
         vC2_2 = vec_madd(vA2, vB2, vC2_2);
         vC3_2 = vec_madd(vA3, vB2, vC3_2); vB2 = vec_ld(0, pB0+KB2+4);
         vC0_3 = vec_madd(vA0, vB3, vC0_3); vA0 = vec_ld(0, pA0+4);
         vC1_3 = vec_madd(vA1, vB3, vC1_3); vA1 = vec_ld(0, pA0+KB+4);
         vC2_3 = vec_madd(vA2, vB3, vC2_3); vA2 = vec_ld(0, pA0+KB2+4);
         vC3_3 = vec_madd(vA3, vB3, vC3_3); vA3 = vec_ld(0, pA0+KB3+4);

         vC0_0 = vec_madd(vA0, vB0, vC0_0); vB3 = vec_ld(0, pB0+KB3+4);
         vC1_0 = vec_madd(vA1, vB0, vC1_0); ATL_pfavR(pA0+KB3+12, cwrdKB, 0);
         vC2_0 = vec_madd(vA2, vB0, vC2_0); ATL_pfavR(pA0+KB4+12, cwrdKB, 1);
         vC3_0 = vec_madd(vA3, vB0, vC3_0); vB0 = vec_ld(0, pB0+8);
         vC0_1 = vec_madd(vA0, vB1, vC0_1); ATL_pfavR(pA0+KB5+12, cwrdKB, 2);
         vC1_1 = vec_madd(vA1, vB1, vC1_1); ATL_pfavR(pA0+KB6+12, cwrdKB, 3);
         vC2_1 = vec_madd(vA2, vB1, vC2_1);
         vC3_1 = vec_madd(vA3, vB1, vC3_1); vB1 = vec_ld(0, pB0+KB+8);
         vC0_2 = vec_madd(vA0, vB2, vC0_2);
         vC1_2 = vec_madd(vA1, vB2, vC1_2);
         vC2_2 = vec_madd(vA2, vB2, vC2_2);
         vC3_2 = vec_madd(vA3, vB2, vC3_2); vB2 = vec_ld(0, pB0+KB2+8);
         vC0_3 = vec_madd(vA0, vB3, vC0_3); vA0 = vec_ld(0, pA0+8);
         vC1_3 = vec_madd(vA1, vB3, vC1_3); vA1 = vec_ld(0, pA0+KB+8);
         vC2_3 = vec_madd(vA2, vB3, vC2_3); vA2 = vec_ld(0, pA0+KB2+8);
         vC3_3 = vec_madd(vA3, vB3, vC3_3); vA3 = vec_ld(0, pA0+KB3+8);

         vC0_0 = vec_madd(vA0, vB0, vC0_0); vB3 = vec_ld(0, pB0+KB3+8);
         vC1_0 = vec_madd(vA1, vB0, vC1_0);
         vC2_0 = vec_madd(vA2, vB0, vC2_0);
         vC3_0 = vec_madd(vA3, vB0, vC3_0);
         vC0_1 = vec_madd(vA0, vB1, vC0_1);
         vC1_1 = vec_madd(vA1, vB1, vC1_1);
         vC2_1 = vec_madd(vA2, vB1, vC2_1);
         vC3_1 = vec_madd(vA3, vB1, vC3_1);
         vC0_2 = vec_madd(vA0, vB2, vC0_2);
         vC1_2 = vec_madd(vA1, vB2, vC1_2);
         vC2_2 = vec_madd(vA2, vB2, vC2_2);
         vC3_2 = vec_madd(vA3, vB2, vC3_2);
         vC0_3 = vec_madd(vA0, vB3, vC0_3);
         vC1_3 = vec_madd(vA1, vB3, vC1_3); pA0 += 12;
         vC2_3 = vec_madd(vA2, vB3, vC2_3);
         vC3_3 = vec_madd(vA3, vB3, vC3_3); pB0 += 12;

         VecReorder(vC0_0, vC1_0, vC2_0, vC3_0);
         VecReorder(vC0_1, vC1_1, vC2_1, vC3_1);
         vC0_0 = vec_add(vC0_0, vC1_0);
         vC0_1 = vec_add(vC0_1, vC1_1);
         vC2_0 = vec_add(vC2_0, vC3_0); ATL_pfavR(pA0+KB3, cwrdKB, 0);
         vC2_1 = vec_add(vC2_1, vC3_1); ATL_pfavR(pA0+KB4, cwrdKB, 1);
         vC0_0 = vec_add(vC0_0, vC2_0); ATL_pfavR(pA0+KB5, cwrdKB, 2);
         vC0_1 = vec_add(vC0_1, vC2_1); ATL_pfavR(pA0+KB6, cwrdKB, 3);

         VecReorder(vC0_2, vC1_2, vC2_2, vC3_2);
         VecReorder(vC0_3, vC1_3, vC2_3, vC3_3);
         vC0_2 = vec_add(vC0_2, vC1_2);
         vC0_3 = vec_add(vC0_3, vC1_3);
         vC2_2 = vec_add(vC2_2, vC3_2);
         vC2_3 = vec_add(vC2_3, vC3_3);
         vC0_2 = vec_add(vC0_2, vC2_2);
         vC0_3 = vec_add(vC0_3, vC2_3);
         vec_st(vC0_0, 0, tC);
         vec_st(vC0_1, 0, tC+4);
         vec_st(vC0_2, 0, tC+8);
         vec_st(vC0_3, 0, tC+12);
         #ifdef BETAX
            rA0 = beta;
            rC0_0 = rC0_0*rA0 + *tC;
            rC1_0 = rC1_0*rA0 + tC[1];
            rC2_0 = rC2_0*rA0 + tC[2];
            rC3_0 = rC3_0*rA0 + tC[3];

            rC0_1 = rC0_1*rA0 + tC[4];
            rC1_1 = rC1_1*rA0 + tC[5];
            rC2_1 = rC2_1*rA0 + tC[6];
            rC3_1 = rC3_1*rA0 + tC[7];

            rC0_2 = rC0_2*rA0 + tC[8];
            rC1_2 = rC1_2*rA0 + tC[9];
            rC2_2 = rC2_2*rA0 + tC[10];
            rC3_2 = rC3_2*rA0 + tC[11];

            rC0_3 = rC0_3*rA0 + tC[12];
            rC1_3 = rC1_3*rA0 + tC[13];
            rC2_3 = rC2_3*rA0 + tC[14];
            rC3_3 = rC3_3*rA0 + tC[15];
         #else
            rC0_0 += *tC; rC1_0 += tC[1]; rC2_0 += tC[2]; rC3_0 += tC[3];
            rC0_1 += tC[4]; rC1_1 += tC[5]; rC2_1 += tC[6]; rC3_1 += tC[7];
            rC0_2 += tC[8]; rC1_2 += tC[9]; rC2_2 += tC[10]; rC3_2 += tC[11];
            rC0_3 += tC[12]; rC1_3 += tC[13]; rC2_3 += tC[14]; rC3_3 += tC[15];
         #endif

         #ifdef TREAL
            *pC0 = rC0_0; pC0[1] = rC1_0; pC0[2] = rC2_0; pC0[3] = rC3_0;
            *pC1 = rC0_1; pC1[1] = rC1_1; pC1[2] = rC2_1; pC1[3] = rC3_1;
            *pC2 = rC0_2; pC2[1] = rC1_2; pC2[2] = rC2_2; pC2[3] = rC3_2;
            *pC3 = rC0_3; pC3[1] = rC1_3; pC3[2] = rC2_3; pC3[3] = rC3_3;
         #else
            *pC0 = rC0_0; pC0[2] = rC1_0; pC0[4] = rC2_0; pC0[6] = rC3_0;
            *pC1 = rC0_1; pC1[2] = rC1_1; pC1[4] = rC2_1; pC1[6] = rC3_1;
            *pC2 = rC0_2; pC2[2] = rC1_2; pC2[4] = rC2_2; pC2[6] = rC3_2;
            *pC3 = rC0_3; pC3[2] = rC1_3; pC3[4] = rC2_3; pC3[6] = rC3_3;
         #endif
         pC0 += incCm;
         pC1 += incCm;
         pC2 += incCm;
         pC3 += incCm;
         pA0 += incAm;
         pB0 += incBm;
      }
      while(pA0 != stM);
      #ifdef ATL_MCLEAN
      }
         if (M != m4)
            ATL_mmcu(M-m4, 4, K, alpha, pA0, lda, pB0, ldb,
                     beta, pC0, ldc, tC, cwrdKB);
      #endif
      pC0 += incCn; pC1 += incCn; pC2 += incCn; pC3 += incCn;
      pA0 += incAn;
      pB0 += incBn;
   }
   while(pB0 != stN);
   #ifdef ATL_MNCLEAN
   }
   else /* did not enter N-loop */
      ATL_mmcu(M, N, K, alpha, A, lda, B, ldb, beta, C, ldc, tC, cwrdKB);
   if (N > 4 && M > 4 && N != n4) /* entered N-loop, must clean N */
      ATL_mmcu(M, N-n4, K, alpha, A, lda, pB0, ldb,
               beta, pC0, ldc, tC, cwrdKB);
   #endif
   free(vC);
}
#ifdef incAm
   #undef incAm
#endif
#ifdef incBn
   #undef incBn
#endif
#ifdef incCm
   #undef incCm
#endif