Changeset 38
- Timestamp:
- 08/28/07 09:37:54 (1 year ago)
- Files:
-
- 2 added
- 5 modified
-
include/fflas-ffpack/fflas.h (modified) (1 diff)
-
include/fflas-ffpack/fflas_ftrmm.inl (modified) (1 diff)
-
include/fflas-ffpack/fflas_ftrmm_src.inl (modified) (6 diffs)
-
include/fflas-ffpack/fflas_ftrsm.inl (modified) (2 diffs)
-
include/fflas-ffpack/fflas_ftrsm_src.inl (modified) (10 diffs)
-
tests/test-ftrsm.C (added)
-
tests/testeur_ftrsm.C (added)
Legend:
- Unmodified
- Added
- Removed
-
include/fflas-ffpack/fflas.h
r36 r38 620 620 class ftrmmRightLowerTransUnit; 621 621 622 // template<class Field>623 // static void ftrsmLeftUpNoTrans (const Field& F, const FFLAS_DIAG Diag,624 // const size_t M, const size_t N,625 // const typename Field::Element alpha,626 // const typename Field::Element * A, const size_t lda,627 // typename Field::Element * B, const size_t ldb);628 629 // template<class Field>630 // static void ftrsmLeftUpTrans (const Field& F, const FFLAS_DIAG Diag,631 // const size_t M, const size_t N,632 // const typename Field::Element alpha,633 // const typename Field::Element * A, const size_t lda,634 // typename Field::Element * B, const size_t ldb);635 636 // template<class Field>637 // static void ftrsmLeftLowNoTrans (const Field& F, const FFLAS_DIAG Diag,638 // const size_t M, const size_t N,639 // const typename Field::Element alpha,640 // typename Field::Element * A, const size_t lda,641 // typename Field::Element * B, const size_t ldb,642 // const size_t nmax);643 // template<class Element>644 // class callFtrsmLeftLowNoTrans;645 646 // template<class Element>647 // class callFtrsmRightUpNoTrans;648 649 // template<class Element>650 // class callFtrmmLeftUpNoTrans;651 652 // template<class Element>653 // class callFtrmmLeftUpTrans;654 655 // template<class Element>656 // class callFtrmmLeftLowNoTrans;657 658 // template<class Element>659 // class callFtrmmLeftLowTrans;660 661 // template<class Element>662 // class callFtrmmRightUpNoTrans;663 664 // template<class Element>665 // class callFtrmmRightUpTrans;666 667 // template<class Element>668 // class callFtrmmRightLowNoTrans;669 670 // template<class Element>671 // class callFtrmmRightLowTrans;672 673 // template<class Field>674 // static void ftrsmLeftLowTrans (const Field& F, const FFLAS_DIAG Diag,675 // const size_t M, const size_t N,676 // const typename Field::Element alpha,677 // const typename Field::Element * A, const size_t lda,678 // typename Field::Element * B, const size_t ldb);679 680 // template<class Field>681 // static void ftrsmRightUpNoTrans (const Field& F, const FFLAS_DIAG Diag,682 // const size_t M, const size_t N,683 // const typename Field::Element alpha,684 // typename Field::Element * A, const size_t lda,685 // typename Field::Element * B, const size_t ldb,686 // const size_t nmax);687 688 // template<class Field>689 // static void ftrsmRightUpTrans (const Field& F, const FFLAS_DIAG Diag,690 // const size_t M, const size_t N,691 // const typename Field::Element alpha,692 // const typename Field::Element * A, const size_t lda,693 // typename Field::Element * B, const size_t ldb);694 695 // template<class Field>696 // static void ftrsmRightLowNoTrans (const Field& F, const FFLAS_DIAG Diag,697 // const size_t M, const size_t N,698 // const typename Field::Element alpha,699 // const typename Field::Element * A, const size_t lda,700 // typename Field::Element * B, const size_t ldb);701 702 // template<class Field>703 // static void ftrsmRightLowTrans (const Field& F, const FFLAS_DIAG Diag,704 // const size_t M, const size_t N,705 // const typename Field::Element alpha,706 // const typename Field::Element * A, const size_t lda,707 // typename Field::Element * B, const size_t ldb);708 709 // Specialized routines for ftrmm710 // template<class Field>711 // static void ftrmmLeftUpNoTrans (const Field& F, const FFLAS_DIAG Diag,712 // const size_t M, const size_t N,713 // const typename Field::Element * A, const size_t lda,714 // typename Field::Element * B, const size_t ldb,715 // const size_t nmax);716 717 // template<class Field>718 // static void ftrmmLeftUpTrans (const Field& F, const FFLAS_DIAG Diag,719 // const size_t M, const size_t N,720 // const typename Field::Element * A, const size_t lda,721 // typename Field::Element * B, const size_t ldb,722 // const size_t nmax);723 724 // template<class Field>725 // static void ftrmmLeftLowNoTrans (const Field& F, const FFLAS_DIAG Diag,726 // const size_t M, const size_t N,727 // const typename Field::Element * A, const size_t lda,728 // typename Field::Element * B, const size_t ldb,729 // const size_t nmax);730 731 // template<class Field>732 // static void ftrmmLeftLowTrans (const Field& F, const FFLAS_DIAG Diag,733 // const size_t M, const size_t N,734 // const typename Field::Element * A, const size_t lda,735 // typename Field::Element * B, const size_t ldb,736 // const size_t nmax);737 738 // template<class Field>739 // static void ftrmmRightUpNoTrans (const Field& F, const FFLAS_DIAG Diag,740 // const size_t M, const size_t N,741 // const typename Field::Element * A, const size_t lda,742 // typename Field::Element * B, const size_t ldb,743 // const size_t nmax);744 745 // template<class Field>746 // static void ftrmmRightUpTrans (const Field& F, const FFLAS_DIAG Diag,747 // const size_t M, const size_t N,748 // const typename Field::Element * A, const size_t lda,749 // typename Field::Element * B, const size_t ldb,750 // const size_t nmax);751 752 // template<class Field>753 // static void ftrmmRightLowNoTrans (const Field& F, const FFLAS_DIAG Diag,754 // const size_t M, const size_t N,755 // const typename Field::Element * A, const size_t lda,756 // typename Field::Element * B, const size_t ldb,757 // const size_t nmax);758 // template<class Field>759 // static void ftrmmRightLowTrans (const Field& F, const FFLAS_DIAG Diag,760 // const size_t M, const size_t N,761 // const typename Field::Element * A, const size_t lda,762 // typename Field::Element * B, const size_t ldb,763 // const size_t nmax);764 622 }; // class FFLAS 765 623 -
include/fflas-ffpack/fflas_ftrmm.inl
r36 r38 14 14 // Computes B <- alpha.op(A).B, B <- alpha.B.op(A) 15 15 // B is M*N, A is M*M if Side==FflasLeft, N*N if Side==FflasRight 16 //--------------------------------------------------------------------- 16 // Warning : unsafe with Trans == FflasTrans (debugging in progress) 17 // //--------------------------------------------------------------------- 17 18 template<class Field> 18 19 inline void -
include/fflas-ffpack/fflas_ftrmm_src.inl
r36 r38 29 29 #define __FFLAS__Aupdate __FFLAS__Atriang + nsplit * __FFLAS__Arowinc 30 30 #define __FFLAS__Arest A + nbblocsplit * nsplit * (lda+1) 31 #define __FFLAS__B updateB + (nbblocsplit - (i+1)) * nsplit * ldb32 #define __FFLAS__B recB + (nbblocsplit - i) * nsplit * ldb31 #define __FFLAS__Brec B + (nbblocsplit - (i+1)) * nsplit * ldb 32 #define __FFLAS__Bupdate B + (nbblocsplit - i) * nsplit * ldb 33 33 #define __FFLAS__Brest B + nbblocsplit * nsplit * ldb 34 #define __FFLAS__A1 A + nsplit* (lda + 1)35 #define __FFLAS__A2 A + nsplit* __FFLAS__Arowinc34 #define __FFLAS__A1 A + (nsplit) * (lda + 1) 35 #define __FFLAS__A2 A + (nsplit) * __FFLAS__Arowinc 36 36 #define __FFLAS__A3 A 37 #define __FFLAS__B1 B + nsplit * ldb37 #define __FFLAS__B1 B + (nsplit) * ldb 38 38 #define __FFLAS__B2 B 39 39 #else 40 #define __FFLAS__Atriang A + i * nsplit* (lda + 1)41 #define __FFLAS__Aupdate A + ( i+1) * nsplit* __FFLAS__Acolinc40 #define __FFLAS__Atriang A + (nrestsplit + i * nsplit) * (lda + 1) 41 #define __FFLAS__Aupdate A + (nrestsplit + i * nsplit) * __FFLAS__Acolinc 42 42 #define __FFLAS__Arest A 43 #define __FFLAS__B updateB + (nrestsplit + i * nsplit) * ldb44 #define __FFLAS__B rec B + MAX(0, int (nrestsplit) + (i-1) * nsplit) * ldb43 #define __FFLAS__Brec B + (nrestsplit + i * nsplit) * ldb 44 #define __FFLAS__Bupdate B 45 45 #define __FFLAS__Brest B 46 46 #define __FFLAS__A1 A 47 #define __FFLAS__A2 A + nsplit* __FFLAS__Acolinc48 #define __FFLAS__A3 A + nsplit* (lda + 1)49 #define __FFLAS__B1 B 50 #define __FFLAS__B2 B + nsplit* ldb47 #define __FFLAS__A2 A + (M-nsplit) * __FFLAS__Acolinc 48 #define __FFLAS__A3 A + (M-nsplit) * (lda + 1) 49 #define __FFLAS__B1 B 50 #define __FFLAS__B2 B + (M-nsplit) * ldb 51 51 #endif 52 52 #else … … 67 67 #define __FFLAS__Aupdate __FFLAS__Atriang + nsplit * __FFLAS__Acolinc 68 68 #define __FFLAS__Arest A + nbblocsplit * nsplit * (lda+1) 69 #define __FFLAS__B updateB + (nbblocsplit - (i+1)) * nsplit70 #define __FFLAS__B recB + (nbblocsplit - i) * nsplit69 #define __FFLAS__Brec B + (nbblocsplit - (i+1)) * nsplit 70 #define __FFLAS__Bupdate B + (nbblocsplit - i) * nsplit 71 71 #define __FFLAS__Brest B + nbblocsplit * nsplit 72 #define __FFLAS__A1 A + nsplit* (lda + 1)73 #define __FFLAS__A2 A + nsplit* __FFLAS__Acolinc72 #define __FFLAS__A1 A + (nsplit) * (lda + 1) 73 #define __FFLAS__A2 A + (nsplit) * __FFLAS__Acolinc 74 74 #define __FFLAS__A3 A 75 75 #define __FFLAS__B1 B + nsplit 76 #define __FFLAS__B2 B 77 #else 78 #define __FFLAS__Atriang A + i * nsplit * (lda + 1)79 #define __FFLAS__Aupdate A + ( i+1) * nsplit * __FFLAS__Arowinc76 #define __FFLAS__B2 B 77 #else 78 #define __FFLAS__Atriang A + (nrestsplit + i * nsplit) * (lda + 1) 79 #define __FFLAS__Aupdate A + (nrestsplit + i * nsplit) * __FFLAS__Arowinc 80 80 #define __FFLAS__Arest A 81 #define __FFLAS__B updateB + (nrestsplit + i * nsplit)82 #define __FFLAS__B rec B + MAX(0, int (nrestsplit) + (i-1) * nsplit)81 #define __FFLAS__Brec B + (nrestsplit + i * nsplit) 82 #define __FFLAS__Bupdate B 83 83 #define __FFLAS__Brest B 84 84 #define __FFLAS__A1 A 85 #define __FFLAS__A2 A + nsplit* __FFLAS__Arowinc86 #define __FFLAS__A3 A + nsplit* (lda + 1)87 #define __FFLAS__B1 B 88 #define __FFLAS__B2 B + nsplit85 #define __FFLAS__A2 A + (N-nsplit) * __FFLAS__Arowinc 86 #define __FFLAS__A3 A + (N-nsplit) * (lda + 1) 87 #define __FFLAS__B1 B 88 #define __FFLAS__B2 B + N-nsplit 89 89 #endif 90 90 #endif … … 159 159 F.neg(Mone, one); 160 160 161 static __FFLAS__DOMAIN D;162 //size_t nblas = TRSMBound<Field> (F);163 164 //size_t nbblocsblas = __FFLAS__Na / nblas;165 //size_t nrestblas = __FFLAS__Na % nblas;166 161 size_t nsplit = DotProdBound (F, 0, one, 167 162 #ifdef __FFLAS__DOUBLE … … 171 166 #endif 172 167 ); 173 //ndel = (ndel / nblas)*nblas; 174 //size_t nsplit = ndel;//MIN (ndel, (nbblocsblas+1) / 2 * nblas); 168 175 169 size_t nbblocsplit = (__FFLAS__Na-1) / nsplit; 176 170 size_t nrestsplit = ((__FFLAS__Na-1) % nsplit) +1; 177 171 178 //std::cout<<"nblas, ndel, nsplit, nbblocsplit, nrestsplit = "<<nblas<<" "<<ndel<<" " 179 //<<nsplit<<" "<<nbblocsplit<<" "<<nrestsplit<<std::endl; 180 181 if (nrestsplit){ 182 //std::cerr<<"nblas nrestsplit, M, N = "<<nblas<<" "<<nrestsplit<<" "<<M<<" "<<N<<std::endl; 172 if (nrestsplit) 183 173 this->delayed (F, __FFLAS__Mbrest, __FFLAS__Nbrest, 184 174 __FFLAS__Arest, lda, __FFLAS__Brest, ldb); 185 } 186 175 187 176 for ( size_t i = 0; i < nbblocsplit; ++i) { 188 177 189 // cerr<<"M,N,K = "<<M<<" "<<(N-(i+1)*nsplit)<<" "<<nsplit<<endl;190 178 #ifdef __FFLAS__RIGHT 191 179 fgemm (F, FflasNoTrans, Mjoin (Fflas, __FFLAS__TRANS), … … 197 185 __FFLAS__Aupdate, lda, __FFLAS__Brec, ldb, one, __FFLAS__Bupdate, ldb); 198 186 #endif 187 199 188 this->delayed (F, __FFLAS__Mb, __FFLAS__Nb, 200 189 __FFLAS__Atriang, lda, __FFLAS__Brec, ldb); … … 222 211 F.neg(Mone,one); 223 212 if (__FFLAS__Na == 1) 213 #ifdef __FFLAS__NONUNIT 224 214 fscal(F, __FFLAS__Bdim, *A, B, __FFLAS__Bnorminc); 225 226 else { // __FFLAS__Na > 1 215 #else 216 ; 217 #endif 218 219 else { // __FFLAS__Na > 1 227 220 size_t nsplit = __FFLAS__Na >> 1; 228 this->operator() (F, __FFLAS__Mb, __FFLAS__Nb, __FFLAS__A1, lda, __FFLAS__B1, ldb); 229 // cerr<<"delay, Nup = "<<delay<<" "<<Nup<<endl; 230 //cerr<<"M,N,K = "<<M<<" "<<Ndown<<" "<<Nup<<endl; 231 221 this->operator() (F, __FFLAS__Mb2, __FFLAS__Nb2, __FFLAS__A1, lda, __FFLAS__B1, ldb); 222 232 223 #ifdef __FFLAS__RIGHT 233 224 fgemm (F, FflasNoTrans , Mjoin (Fflas, __FFLAS__TRANS), 234 225 __FFLAS__Mb2, __FFLAS__Nb2, nsplit, one, 235 __FFLAS__B 1, ldb, __FFLAS__A2, lda, one, __FFLAS__B2, ldb);226 __FFLAS__B2, ldb, __FFLAS__A2, lda, one, __FFLAS__B1, ldb); 236 227 #else 237 228 fgemm (F, Mjoin (Fflas, __FFLAS__TRANS), FflasNoTrans, 238 229 __FFLAS__Mb2, __FFLAS__Nb2, nsplit, one, 239 __FFLAS__A2, lda, __FFLAS__B 1, ldb, one, __FFLAS__B2, ldb);240 #endif 241 this->operator() (F, __FFLAS__Mb 2, __FFLAS__Nb2, __FFLAS__A3, lda, __FFLAS__B2, ldb);230 __FFLAS__A2, lda, __FFLAS__B2, ldb, one, __FFLAS__B1, ldb); 231 #endif 232 this->operator() (F, __FFLAS__Mb, __FFLAS__Nb, __FFLAS__A3, lda, __FFLAS__B2, ldb); 242 233 } 243 234 } -
include/fflas-ffpack/fflas_ftrsm.inl
r36 r38 14 14 // Computes B <- alpha.op(A^-1).B, B <- alpha.B.op(A^-1) 15 15 // B is M*N, A is M*M if Side==FflasLeft, N*N if Side==FflasRight 16 // Warning : unsafe with Trans == FflasTrans (debugging in progress) 16 17 //--------------------------------------------------------------------- 17 18 template<class Field> … … 678 679 #undef __FFLAS__TRANSPOSE 679 680 #undef __FFLAS__UNIT 680 681 //682 // template<class Field>683 // inline void684 // FFLAS::ftrsmLeftUpNoTrans(const Field& F, const FFLAS_DIAG Diag,685 // const size_t M, const size_t N,686 // const typename Field::Element alpha,687 // const typename Field::Element * A, const size_t lda,688 // typename Field::Element * B, const size_t ldb){689 // static typename Field::Element Mone;690 // static typename Field::Element one;691 // F.init(Mone, -1.0);692 // F.init(one, 1.0);693 // if ( M==1 ){694 // if (Diag == FflasNonUnit ){695 // typename Field::Element inv;696 // F.inv(inv, *A);697 // fscal(F, N, inv, B, 1);698 // }699 // }700 // else{701 // size_t Mup=M>>1;702 // size_t Mdown = M-Mup;703 // ftrsmLeftUpNoTrans( F, Diag, Mdown, N, alpha,704 // A+Mup*(lda+1), lda, B+Mup*ldb, ldb);705 // fgemm( F, FflasNoTrans, FflasNoTrans, Mup, N, Mdown, Mone,706 // A+Mup, lda, B+Mup*ldb, ldb, alpha, B, ldb);707 // ftrsmLeftUpNoTrans( F, Diag, Mup, N, one, A, lda, B, ldb);708 // }709 710 // }711 712 // template<class Field>713 // inline void714 // FFLAS::ftrsmLeftUpTrans(const Field& F, const FFLAS_DIAG Diag,715 // const size_t M, const size_t N,716 // const typename Field::Element alpha,717 // const typename Field::Element * A, const size_t lda,718 // typename Field::Element * B, const size_t ldb){719 720 // static typename Field::Element Mone;721 // static typename Field::Element one;722 // F.init(Mone, -1.0);723 // F.init(one, 1.0);724 // if ( M==1 ){725 // if (Diag == FflasNonUnit ){726 // typename Field::Element inv;727 // F.inv(inv, *A);728 // fscal(F, N, inv, B, 1);729 // }730 // }731 // else{732 // size_t Mup=M>>1;733 // size_t Mdown = M-Mup;734 // ftrsmLeftUpTrans( F, Diag, Mdown, N, alpha,735 // A+Mup*(lda+1), lda, B+Mup*ldb, ldb);736 // fgemm( F, FflasTrans, FflasNoTrans, Mup, N, Mdown, Mone,737 // A+Mup*lda, lda, B+Mup*ldb, ldb, alpha, B, ldb);738 // ftrsmLeftUpTrans( F, Diag, Mup, N, one, A, lda, B, ldb);739 // }740 // }741 742 // template<class Field>743 // inline void744 // FFLAS::ftrsmLeftLowNoTrans(const Field& F, const FFLAS_DIAG Diag,745 // const size_t M, const size_t N,746 // const typename Field::Element alpha,747 // typename Field::Element * A, const size_t lda,748 // typename Field::Element * B, const size_t ldb, const size_t nmax){749 750 // callFtrsmLeftLowNoTrans<typename Field::Element>() (F,Diag,M,N,alpha,A,lda,B,ldb,nmax);751 // }752 753 754 755 // // Implementation of Ftrsmleftlownotrans on a Field with a756 // // double floating point representation757 // template<>758 // class FFLAS::callFtrsmLeftLowNoTrans<double>{759 // public:760 // template <class Field>761 // void operator ()(const Field& F, const FFLAS_DIAG Diag,762 // const size_t M, const size_t N,763 // const typename Field::Element alpha,764 // typename Field::Element * A, const size_t lda,765 // typename Field::Element * B, const size_t ldb, const size_t nmax){766 767 // static typename Field::Element Mone;768 // static typename Field::Element one;769 // F.init(one, 1.0);770 // F.neg(Mone, one);771 772 // if ( M <= nmax ){773 // typename Field::Element inv;774 // if (Diag == FflasNonUnit ){775 // //Normalization of A and correction of B776 // typename Field::Element * Ai = A;777 // typename Field::Element * Bi = B;778 // for (size_t i=0; i<M; ++i){779 // F.inv( inv, *(Ai+i) );780 // fscal(F, i, inv, Ai, 1 );781 // fscal(F, N, inv, Bi, 1 );782 // Ai += lda; Bi+=ldb;783 784 // }785 // }786 787 // cblas_dtrsm( CblasRowMajor, CblasLeft, CblasLower, CblasNoTrans,788 // CblasUnit, M, N, one, A, lda, B, ldb );789 // for (size_t i=0; i< M; ++i)790 // for (size_t j=0; j<N; ++j)791 // F.init(*(B+i*ldb+j),*(B+i*ldb+j));792 793 794 // if (Diag == FflasNonUnit ){795 // //Denormalization of A796 // typename Field::Element * Ai=A;797 // for (size_t i=0; i<M; ++i){798 // fscal( F, i, *(Ai+i), Ai, 1 );799 // Ai += lda;800 // }801 // }802 // }803 // else{804 // size_t Mup=M>>1;805 // size_t Mdown = M-Mup;806 // this->operator()( F, Diag, Mup, N, alpha, A, lda, B, ldb, nmax);807 // fgemm( F, FflasNoTrans, FflasNoTrans, Mdown, N, Mup,808 // Mone, A+Mup*lda, lda, B, ldb, alpha, B+Mup*ldb, ldb);809 // this->operator()( F, Diag, Mdown, N, one,810 // A+Mup*(lda+1), lda, B+Mup*ldb, ldb, nmax);811 // }812 813 // }814 // };815 816 // template<>817 // class FFLAS::callFtrsmLeftLowNoTrans<float>{818 // public:819 // template <class Field>820 // void operator ()(const Field& F, const FFLAS_DIAG Diag,821 // const size_t M, const size_t N,822 // const typename Field::Element alpha,823 // typename Field::Element * A, const size_t lda,824 // typename Field::Element * B, const size_t ldb, const size_t nmax){825 826 // static typename Field::Element Mone;827 // static typename Field::Element one;828 // F.init(one, 1.0);829 // F.neg(Mone, one);830 831 // if ( M <= nmax ){832 // typename Field::Element inv;833 // if (Diag == FflasNonUnit ){834 // //Normalization of A and correction of B835 // typename Field::Element * Ai = A;836 // typename Field::Element * Bi = B;837 // for (size_t i=0; i<M; ++i){838 // F.inv( inv, *(Ai+i) );839 // fscal(F, i, inv, Ai, 1 );840 // fscal(F, N, inv, Bi, 1 );841 // Ai += lda; Bi+=ldb;842 843 // }844 // }845 846 // cblas_strsm( CblasRowMajor, CblasLeft, CblasLower, CblasNoTrans,847 // CblasUnit, M, N, alpha, A, lda, B, ldb );848 // for (size_t i=0; i< M; ++i)849 // for (size_t j=0; j<N; ++j)850 // F.init(*(B+i*ldb+j),*(B+i*ldb+j));851 852 // if (Diag == FflasNonUnit ){853 // //Denormalization of A854 // typename Field::Element * Ai=A;855 // for (size_t i=0; i<M; ++i){856 // fscal( F, i, *(Ai+i), Ai, 1 );857 // Ai += lda;858 // }859 // }860 // }861 // else{862 // size_t Mup=M>>1;863 // size_t Mdown = M-Mup;864 // this->operator()( F, Diag, Mup, N, alpha, A, lda, B, ldb, nmax);865 // fgemm( F, FflasNoTrans, FflasNoTrans, Mdown, N, Mup,866 // Mone, A+Mup*lda, lda, B, ldb, alpha, B+Mup*ldb, ldb);867 // this->operator()( F, Diag, Mdown, N, one,868 // A+Mup*(lda+1), lda, B+Mup*ldb, ldb, nmax);869 // }870 871 // }872 // };873 874 // template<class Element>875 // class FFLAS::callFtrsmLeftLowNoTrans{876 // public:877 // template <class Field>878 // void operator()(const Field& F, const FFLAS_DIAG Diag,879 // const size_t M, const size_t N,880 // const typename Field::Element alpha,881 // typename Field::Element * A, const size_t lda,882 // typename Field::Element * B, const size_t ldb, const size_t nmax){883 884 // static typename Field::Element Mone;885 // static typename Field::Element one;886 // F.init(Mone, -1.0);887 // F.init(one, 1.0);888 // if ( M <= nmax ){889 // typename Field::Element inv;890 // if (Diag == FflasNonUnit ){891 // //Normalization of A and correction of B892 // // A<-DA, B<-DB893 // typename Field::Element * Ai = A;894 // typename Field::Element * Bi = B;895 // for (size_t i=0; i<M; ++i){896 // F.inv( inv, *(Ai+i) );897 // fscal(F, i, inv, Ai, 1 );898 // fscal(F, N, inv, Bi, 1 );899 // Ai += lda; Bi+=ldb;900 // }901 // }902 // double alphad;903 // if (F.areEqual(alpha, Mone))904 // alphad = -1.0;905 // else906 // F.convert( alphad, alpha );907 // DoubleDomain::Element * Ad = new DoubleDomain::Element[M*M];908 // DoubleDomain::Element * Bd = new DoubleDomain::Element[M*N];909 // MatF2MatD( F, Ad, N, A, lda, M, M );910 // MatF2MatD( F, Bd, N, B, ldb, M, N );911 // cblas_dtrsm( CblasRowMajor, CblasLeft, CblasLower, CblasNoTrans,912 // CblasUnit, M, N, alphad, Ad, M, Bd, N );913 // delete[] Ad;914 // MatD2MatF( F, B, ldb, Bd, N, M, N );915 // delete[] Bd;916 // if (Diag == FflasNonUnit ){917 // //Denormalization of A918 // // A-> D^(-1)A919 // typename Field::Element * Ai=A;920 // for (size_t i=0; i<M; ++i){921 // fscal( F, i, *(Ai+i), Ai, 1 );