diff options
Diffstat (limited to 'lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c')
-rw-r--r-- | lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c | 78 |
1 files changed, 0 insertions, 78 deletions
diff --git a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c index 61cb002..310f6ee 100644 --- a/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c +++ b/lib/Numeric/LinearAlgebra/LAPACK/lapack-aux.c | |||
@@ -814,81 +814,3 @@ int lu_l_C(KCMAT(a), DVEC(ipiv), CMAT(r)) { | |||
814 | free(auxipiv); | 814 | free(auxipiv); |
815 | OK | 815 | OK |
816 | } | 816 | } |
817 | |||
818 | //////////////////////////////////////////////////////////// | ||
819 | |||
820 | // int multiplyR(KDMAT(a),KDMAT(b),DMAT(r)) { | ||
821 | // REQUIRES(ac==br && ar==rr && bc==rc,BAD_SIZE); | ||
822 | // int i,j,k; | ||
823 | // for (i=0;i<ar;i++) { | ||
824 | // for(j=0;j<bc;j++) { | ||
825 | // double temp = 0; | ||
826 | // for(k=0;k<ac;k++) { | ||
827 | // temp += ap[i*ac+k]*bp[k*bc+j]; | ||
828 | // } | ||
829 | // rp[i*rc+j] = temp; | ||
830 | // } | ||
831 | // } | ||
832 | // OK | ||
833 | // } | ||
834 | |||
835 | |||
836 | int multiplyC(KCMAT(a),KCMAT(b),CMAT(r)) { | ||
837 | REQUIRES(ac==br && ar==rr && bc==rc,BAD_SIZE); | ||
838 | int i,j,k; | ||
839 | for (i=0;i<ar;i++) { | ||
840 | for(j=0;j<bc;j++) { | ||
841 | doublecomplex temp = {0,0}; | ||
842 | for(k=0;k<ac;k++) { | ||
843 | doublecomplex A = ((doublecomplex*)ap)[i*ac+k]; | ||
844 | doublecomplex B = ((doublecomplex*)bp)[k*bc+j]; | ||
845 | double w = A.r * B.r - A.i * B.i; | ||
846 | double w1 = A.r * B.r; | ||
847 | double w2 = A.i * B.i; | ||
848 | if(w != w) { | ||
849 | printf("at : %d %d %d\n",i,j,k); | ||
850 | printf("%f %f %f\n",A.i,B.i, A.i * B.i); | ||
851 | printf("%f %f %f\n",A.i,B.i, w2); | ||
852 | exit(1); | ||
853 | } | ||
854 | temp.r += (w + w1-w2)/2; | ||
855 | //temp.r += w; | ||
856 | temp.i += A.r * B.i + A.i * B.r; | ||
857 | } | ||
858 | ((doublecomplex*)rp)[i*rc+j] = temp; | ||
859 | } | ||
860 | } | ||
861 | OK | ||
862 | } | ||
863 | |||
864 | // void dgemm_(char *, char *, integer *, integer *, integer *, | ||
865 | // double *, const double *, integer *, const double *, | ||
866 | // integer *, double *, double *, integer *); | ||
867 | // | ||
868 | // void zgemm_(char *, char *, integer *, integer *, integer *, | ||
869 | // doublecomplex *, const doublecomplex *, integer *, const doublecomplex *, | ||
870 | // integer *, doublecomplex *, doublecomplex *, integer *); | ||
871 | // | ||
872 | // | ||
873 | // int multiplyR2(KDMAT(a),KDMAT(b),DMAT(r)) { | ||
874 | // REQUIRES(ac==br && ar==rr && bc==rc,BAD_SIZE); | ||
875 | // double alpha = 1; | ||
876 | // double beta = 0; | ||
877 | // integer m = ar; | ||
878 | // integer n = bc; | ||
879 | // integer k = ac; | ||
880 | // int i,j; | ||
881 | // dgemm_("N","N",&m,&n,&k,&alpha,ap,&m,bp,&k,&beta,rp,&m); | ||
882 | // OK | ||
883 | // } | ||
884 | // | ||
885 | // int multiplyC2(KCMAT(a),KCMAT(b),CMAT(r)) { | ||
886 | // REQUIRES(ac==br && ar==rr && bc==rc,BAD_SIZE); | ||
887 | // integer m = ar; | ||
888 | // integer n = bc; | ||
889 | // integer k = ac; | ||
890 | // doublecomplex alpha = {1,0}; | ||
891 | // doublecomplex beta = {0,0}; | ||
892 | // zgemm_("N","N",&m,&n,&k,&alpha,(doublecomplex*)ap,&m,(doublecomplex*)bp,&k,&beta,(doublecomplex*)rp,&m); | ||
893 | // OK | ||
894 | // } | ||