| 1041 | } |
| 1042 | |
| 1043 | void BaseVectorPtrMV :: Add (const MultiVector & v2, FlatMatrix<double> mat) |
| 1044 | { |
| 1045 | static Timer t("BaseVector-MV :: mult mat"); |
| 1046 | RegionTimer reg(t); |
| 1047 | t.AddFlops (mat.Height()*mat.Width()*this->RefVec()->FVDouble().Size()); |
| 1048 | |
| 1049 | size_t n = refvec->FVDouble().Size(); |
| 1050 | |
| 1051 | size_t BBH = 256; |
| 1052 | size_t AH = 512; |
| 1053 | size_t BH = 128; |
| 1054 | |
| 1055 | ParallelFor(1 + n / BBH, [&] (int i) { |
| 1056 | |
| 1057 | int i0 = BBH * i; |
| 1058 | int is = min(BBH, n - i0); |
| 1059 | |
| 1060 | // allocate memory for pointer arrays |
| 1061 | STACK_ARRAY(double*, ppx, AH); |
| 1062 | STACK_ARRAY(double*, ppy, BH); |
| 1063 | |
| 1064 | for (int j0 = 0; j0 < Size(); j0 += AH) { |
| 1065 | int js = min(AH, Size() - j0); |
| 1066 | |
| 1067 | // get pointers of first multivector |
| 1068 | for (int ell = 0; ell < js; ell++) { |
| 1069 | ppx[ell] = (*this)[j0 + ell]->FVDouble().Addr(i0); |
| 1070 | } |
| 1071 | |
| 1072 | for (int k0 = 0; k0 < v2.Size(); k0+=BH) { |
| 1073 | int ks = min(BH, v2.Size() - k0); |
| 1074 | |
| 1075 | // get pointers of second multivector |
| 1076 | for (int ell = 0; ell < ks; ell++) { |
| 1077 | ppy[ell] = v2[k0 + ell]->FVDouble().Addr(i0); |
| 1078 | } |
| 1079 | |
| 1080 | MultiVectorAdd(is, FlatArray(js, ppx), FlatArray(ks, ppy), SliceMatrix<double>(ks, js, mat.Width(), &mat(k0,j0))); |
| 1081 | } |
| 1082 | |
| 1083 | } |
| 1084 | |
| 1085 | }); |
| 1086 | |
| 1087 | } |
| 1088 | |
| 1089 | |
| 1090 | void BaseVectorPtrMV :: Add (const MultiVector & v2, FlatMatrix<Complex> mat) |