Mercurial > hg > octave-lyh

comparison liboctave/CSparse.cc @ 5785:6b9cec830d72

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 2006-05-03 19:32:46 by dbateman]
author dbateman
date Wed, 03 May 2006 19:32:48 +0000
parents cbf717bf8150
children 565d0cd4d9d0
comparison
equal deleted inserted replaced
5784:70f7659d0fb9 5785:6b9cec830d72
605 SparseComplexMatrix 605 SparseComplexMatrix
606 SparseComplexMatrix::inverse (void) const 606 SparseComplexMatrix::inverse (void) const
607 { 607 {
608 octave_idx_type info; 608 octave_idx_type info;
609 double rcond; 609 double rcond;
610 SparseType mattype (*this); 610 MatrixType mattype (*this);
611 return inverse (mattype, info, rcond, 0, 0); 611 return inverse (mattype, info, rcond, 0, 0);
612 } 612 }
613 613
614 SparseComplexMatrix 614 SparseComplexMatrix
615 SparseComplexMatrix::inverse (SparseType& mattype) const 615 SparseComplexMatrix::inverse (MatrixType& mattype) const
616 { 616 {
617 octave_idx_type info; 617 octave_idx_type info;
618 double rcond; 618 double rcond;
619 return inverse (mattype, info, rcond, 0, 0); 619 return inverse (mattype, info, rcond, 0, 0);
620 } 620 }
621 621
622 SparseComplexMatrix 622 SparseComplexMatrix
623 SparseComplexMatrix::inverse (SparseType& mattype, octave_idx_type& info) const 623 SparseComplexMatrix::inverse (MatrixType& mattype, octave_idx_type& info) const
624 { 624 {
625 double rcond; 625 double rcond;
626 return inverse (mattype, info, rcond, 0, 0); 626 return inverse (mattype, info, rcond, 0, 0);
627 } 627 }
628 628
629 SparseComplexMatrix 629 SparseComplexMatrix
630 SparseComplexMatrix::dinverse (SparseType &mattyp, octave_idx_type& info, 630 SparseComplexMatrix::dinverse (MatrixType &mattyp, octave_idx_type& info,
631 double& rcond, const bool, 631 double& rcond, const bool,
632 const bool calccond) const 632 const bool calccond) const
633 { 633 {
634 SparseComplexMatrix retval; 634 SparseComplexMatrix retval;
635 635
643 { 643 {
644 // Print spparms("spumoni") info if requested 644 // Print spparms("spumoni") info if requested
645 int typ = mattyp.type (); 645 int typ = mattyp.type ();
646 mattyp.info (); 646 mattyp.info ();
647 647
648 if (typ == SparseType::Diagonal || 648 if (typ == MatrixType::Diagonal ||
649 typ == SparseType::Permuted_Diagonal) 649 typ == MatrixType::Permuted_Diagonal)
650 { 650 {
651 if (typ == SparseType::Permuted_Diagonal) 651 if (typ == MatrixType::Permuted_Diagonal)
652 retval = transpose(); 652 retval = transpose();
653 else 653 else
654 retval = *this; 654 retval = *this;
655 655
656 // Force make_unique to be called 656 // Force make_unique to be called
679 679
680 return retval; 680 return retval;
681 } 681 }
682 682
683 SparseComplexMatrix 683 SparseComplexMatrix
684 SparseComplexMatrix::tinverse (SparseType &mattyp, octave_idx_type& info, 684 SparseComplexMatrix::tinverse (MatrixType &mattyp, octave_idx_type& info,
685 double& rcond, const bool, 685 double& rcond, const bool,
686 const bool calccond) const 686 const bool calccond) const
687 { 687 {
688 SparseComplexMatrix retval; 688 SparseComplexMatrix retval;
689 689
697 { 697 {
698 // Print spparms("spumoni") info if requested 698 // Print spparms("spumoni") info if requested
699 int typ = mattyp.type (); 699 int typ = mattyp.type ();
700 mattyp.info (); 700 mattyp.info ();
701 701
702 if (typ == SparseType::Upper || typ == SparseType::Permuted_Upper || 702 if (typ == MatrixType::Upper || typ == MatrixType::Permuted_Upper ||
703 typ == SparseType::Lower || typ == SparseType::Permuted_Lower) 703 typ == MatrixType::Lower || typ == MatrixType::Permuted_Lower)
704 { 704 {
705 double anorm = 0.; 705 double anorm = 0.;
706 double ainvnorm = 0.; 706 double ainvnorm = 0.;
707 707
708 if (calccond) 708 if (calccond)
716 if (atmp > anorm) 716 if (atmp > anorm)
717 anorm = atmp; 717 anorm = atmp;
718 } 718 }
719 } 719 }
720 720
721 if (typ == SparseType::Upper || typ == SparseType::Lower) 721 if (typ == MatrixType::Upper || typ == MatrixType::Lower)
722 { 722 {
723 octave_idx_type nz = nnz (); 723 octave_idx_type nz = nnz ();
724 octave_idx_type cx = 0; 724 octave_idx_type cx = 0;
725 octave_idx_type nz2 = nz; 725 octave_idx_type nz2 = nz;
726 retval = SparseComplexMatrix (nr, nc, nz2); 726 retval = SparseComplexMatrix (nr, nc, nz2);
812 812
813 OCTAVE_LOCAL_BUFFER (Complex, work, nr); 813 OCTAVE_LOCAL_BUFFER (Complex, work, nr);
814 OCTAVE_LOCAL_BUFFER (octave_idx_type, rperm, nr); 814 OCTAVE_LOCAL_BUFFER (octave_idx_type, rperm, nr);
815 815
816 octave_idx_type *perm = mattyp.triangular_perm(); 816 octave_idx_type *perm = mattyp.triangular_perm();
817 if (typ == SparseType::Permuted_Upper) 817 if (typ == MatrixType::Permuted_Upper)
818 { 818 {
819 for (octave_idx_type i = 0; i < nr; i++) 819 for (octave_idx_type i = 0; i < nr; i++)
820 rperm[perm[i]] = i; 820 rperm[perm[i]] = i;
821 } 821 }
822 else 822 else
917 917
918 return retval; 918 return retval;
919 } 919 }
920 920
921 SparseComplexMatrix 921 SparseComplexMatrix
922 SparseComplexMatrix::inverse (SparseType& mattype, octave_idx_type& info, 922 SparseComplexMatrix::inverse (MatrixType& mattype, octave_idx_type& info,
923 double& rcond, int, int calc_cond) const 923 double& rcond, int, int calc_cond) const
924 { 924 {
925 int typ = mattype.type (false); 925 int typ = mattype.type (false);
926 SparseComplexMatrix ret; 926 SparseComplexMatrix ret;
927 927
928 if (typ == SparseType::Unknown) 928 if (typ == MatrixType::Unknown)
929 typ = mattype.type (*this); 929 typ = mattype.type (*this);
930 930
931 if (typ == SparseType::Diagonal || typ == SparseType::Permuted_Diagonal) 931 if (typ == MatrixType::Diagonal || typ == MatrixType::Permuted_Diagonal)
932 ret = dinverse (mattype, info, rcond, true, calc_cond); 932 ret = dinverse (mattype, info, rcond, true, calc_cond);
933 else if (typ == SparseType::Upper || typ == SparseType::Permuted_Upper) 933 else if (typ == MatrixType::Upper || typ == MatrixType::Permuted_Upper)
934 ret = tinverse (mattype, info, rcond, true, calc_cond).transpose(); 934 ret = tinverse (mattype, info, rcond, true, calc_cond).transpose();
935 else if (typ == SparseType::Lower || typ == SparseType::Permuted_Lower) 935 else if (typ == MatrixType::Lower || typ == MatrixType::Permuted_Lower)
936 ret = transpose().tinverse (mattype, info, rcond, true, calc_cond); 936 ret = transpose().tinverse (mattype, info, rcond, true, calc_cond);
937 else if (typ != SparseType::Rectangular) 937 else if (typ != MatrixType::Rectangular)
938 { 938 {
939 if (mattype.is_hermitian()) 939 if (mattype.is_hermitian())
940 { 940 {
941 SparseType tmp_typ (SparseType::Upper); 941 MatrixType tmp_typ (MatrixType::Upper);
942 SparseComplexCHOL fact (*this, info, false); 942 SparseComplexCHOL fact (*this, info, false);
943 rcond = fact.rcond(); 943 rcond = fact.rcond();
944 if (info == 0) 944 if (info == 0)
945 { 945 {
946 double rcond2; 946 double rcond2;
951 } 951 }
952 else 952 else
953 { 953 {
954 // Matrix is either singular or not positive definite 954 // Matrix is either singular or not positive definite
955 mattype.mark_as_unsymmetric (); 955 mattype.mark_as_unsymmetric ();
956 typ = SparseType::Full; 956 typ = MatrixType::Full;
957 } 957 }
958 } 958 }
959 959
960 if (!mattype.is_hermitian()) 960 if (!mattype.is_hermitian())
961 { 961 {
962 octave_idx_type n = rows(); 962 octave_idx_type n = rows();
963 ColumnVector Qinit(n); 963 ColumnVector Qinit(n);
964 for (octave_idx_type i = 0; i < n; i++) 964 for (octave_idx_type i = 0; i < n; i++)
965 Qinit(i) = i; 965 Qinit(i) = i;
966 966
967 SparseType tmp_typ (SparseType::Upper); 967 MatrixType tmp_typ (MatrixType::Upper);
968 SparseComplexLU fact (*this, Qinit, -1.0, false); 968 SparseComplexLU fact (*this, Qinit, -1.0, false);
969 rcond = fact.rcond(); 969 rcond = fact.rcond();
970 double rcond2; 970 double rcond2;
971 SparseComplexMatrix InvL = fact.L().transpose(). 971 SparseComplexMatrix InvL = fact.L().transpose().
972 tinverse(tmp_typ, info, rcond2, true, false); 972 tinverse(tmp_typ, info, rcond2, true, false);
1123 1123
1124 return retval; 1124 return retval;
1125 } 1125 }
1126 1126
1127 ComplexMatrix 1127 ComplexMatrix
1128 SparseComplexMatrix::dsolve (SparseType &mattype, const Matrix& b, 1128 SparseComplexMatrix::dsolve (MatrixType &mattype, const Matrix& b,
1129 octave_idx_type& err, double& rcond, 1129 octave_idx_type& err, double& rcond,
1130 solve_singularity_handler, bool calc_cond) const 1130 solve_singularity_handler, bool calc_cond) const
1131 { 1131 {
1132 ComplexMatrix retval; 1132 ComplexMatrix retval;
1133 1133
1143 { 1143 {
1144 // Print spparms("spumoni") info if requested 1144 // Print spparms("spumoni") info if requested
1145 int typ = mattype.type (); 1145 int typ = mattype.type ();
1146 mattype.info (); 1146 mattype.info ();
1147 1147
1148 if (typ == SparseType::Diagonal || 1148 if (typ == MatrixType::Diagonal ||
1149 typ == SparseType::Permuted_Diagonal) 1149 typ == MatrixType::Permuted_Diagonal)
1150 { 1150 {
1151 retval.resize (nc, b.cols(), Complex(0.,0.)); 1151 retval.resize (nc, b.cols(), Complex(0.,0.));
1152 if (typ == SparseType::Diagonal) 1152 if (typ == MatrixType::Diagonal)
1153 for (octave_idx_type j = 0; j < b.cols(); j++) 1153 for (octave_idx_type j = 0; j < b.cols(); j++)
1154 for (octave_idx_type i = 0; i < nm; i++) 1154 for (octave_idx_type i = 0; i < nm; i++)
1155 retval(i,j) = b(i,j) / data (i); 1155 retval(i,j) = b(i,j) / data (i);
1156 else 1156 else
1157 for (octave_idx_type j = 0; j < b.cols(); j++) 1157 for (octave_idx_type j = 0; j < b.cols(); j++)
1181 1181
1182 return retval; 1182 return retval;
1183 } 1183 }
1184 1184
1185 SparseComplexMatrix 1185 SparseComplexMatrix
1186 SparseComplexMatrix::dsolve (SparseType &mattype, const SparseMatrix& b, 1186 SparseComplexMatrix::dsolve (MatrixType &mattype, const SparseMatrix& b,
1187 octave_idx_type& err, double& rcond, 1187 octave_idx_type& err, double& rcond,
1188 solve_singularity_handler, 1188 solve_singularity_handler,
1189 bool calc_cond) const 1189 bool calc_cond) const
1190 { 1190 {
1191 SparseComplexMatrix retval; 1191 SparseComplexMatrix retval;
1202 { 1202 {
1203 // Print spparms("spumoni") info if requested 1203 // Print spparms("spumoni") info if requested
1204 int typ = mattype.type (); 1204 int typ = mattype.type ();
1205 mattype.info (); 1205 mattype.info ();
1206 1206
1207 if (typ == SparseType::Diagonal || 1207 if (typ == MatrixType::Diagonal ||
1208 typ == SparseType::Permuted_Diagonal) 1208 typ == MatrixType::Permuted_Diagonal)
1209 { 1209 {
1210 octave_idx_type b_nc = b.cols (); 1210 octave_idx_type b_nc = b.cols ();
1211 octave_idx_type b_nz = b.nnz (); 1211 octave_idx_type b_nz = b.nnz ();
1212 retval = SparseComplexMatrix (nc, b_nc, b_nz); 1212 retval = SparseComplexMatrix (nc, b_nc, b_nz);
1213 1213
1214 retval.xcidx(0) = 0; 1214 retval.xcidx(0) = 0;
1215 octave_idx_type ii = 0; 1215 octave_idx_type ii = 0;
1216 if (typ == SparseType::Diagonal) 1216 if (typ == MatrixType::Diagonal)
1217 for (octave_idx_type j = 0; j < b.cols(); j++) 1217 for (octave_idx_type j = 0; j < b.cols(); j++)
1218 { 1218 {
1219 for (octave_idx_type i = b.cidx(j); i < b.cidx(j+1); i++) 1219 for (octave_idx_type i = b.cidx(j); i < b.cidx(j+1); i++)
1220 { 1220 {
1221 if (b.ridx(i) >= nm) 1221 if (b.ridx(i) >= nm)
1270 1270
1271 return retval; 1271 return retval;
1272 } 1272 }
1273 1273
1274 ComplexMatrix 1274 ComplexMatrix
1275 SparseComplexMatrix::dsolve (SparseType &mattype, const ComplexMatrix& b, 1275 SparseComplexMatrix::dsolve (MatrixType &mattype, const ComplexMatrix& b,
1276 octave_idx_type& err, double& rcond, 1276 octave_idx_type& err, double& rcond,
1277 solve_singularity_handler, 1277 solve_singularity_handler,
1278 bool calc_cond) const 1278 bool calc_cond) const
1279 { 1279 {
1280 ComplexMatrix retval; 1280 ComplexMatrix retval;
1291 { 1291 {
1292 // Print spparms("spumoni") info if requested 1292 // Print spparms("spumoni") info if requested
1293 int typ = mattype.type (); 1293 int typ = mattype.type ();
1294 mattype.info (); 1294 mattype.info ();
1295 1295
1296 if (typ == SparseType::Diagonal || 1296 if (typ == MatrixType::Diagonal ||
1297 typ == SparseType::Permuted_Diagonal) 1297 typ == MatrixType::Permuted_Diagonal)
1298 { 1298 {
1299 retval.resize (nc, b.cols(), Complex(0.,0.)); 1299 retval.resize (nc, b.cols(), Complex(0.,0.));
1300 if (typ == SparseType::Diagonal) 1300 if (typ == MatrixType::Diagonal)
1301 for (octave_idx_type j = 0; j < b.cols(); j++) 1301 for (octave_idx_type j = 0; j < b.cols(); j++)
1302 for (octave_idx_type i = 0; i < nm; i++) 1302 for (octave_idx_type i = 0; i < nm; i++)
1303 retval(i,j) = b(i,j) / data (i); 1303 retval(i,j) = b(i,j) / data (i);
1304 else 1304 else
1305 for (octave_idx_type j = 0; j < b.cols(); j++) 1305 for (octave_idx_type j = 0; j < b.cols(); j++)
1329 1329
1330 return retval; 1330 return retval;
1331 } 1331 }
1332 1332
1333 SparseComplexMatrix 1333 SparseComplexMatrix
1334 SparseComplexMatrix::dsolve (SparseType &mattype, const SparseComplexMatrix& b, 1334 SparseComplexMatrix::dsolve (MatrixType &mattype, const SparseComplexMatrix& b,
1335 octave_idx_type& err, double& rcond, 1335 octave_idx_type& err, double& rcond,
1336 solve_singularity_handler, 1336 solve_singularity_handler,
1337 bool calc_cond) const 1337 bool calc_cond) const
1338 { 1338 {
1339 SparseComplexMatrix retval; 1339 SparseComplexMatrix retval;
1350 { 1350 {
1351 // Print spparms("spumoni") info if requested 1351 // Print spparms("spumoni") info if requested
1352 int typ = mattype.type (); 1352 int typ = mattype.type ();
1353 mattype.info (); 1353 mattype.info ();
1354 1354
1355 if (typ == SparseType::Diagonal || 1355 if (typ == MatrixType::Diagonal ||
1356 typ == SparseType::Permuted_Diagonal) 1356 typ == MatrixType::Permuted_Diagonal)
1357 { 1357 {
1358 octave_idx_type b_nc = b.cols (); 1358 octave_idx_type b_nc = b.cols ();
1359 octave_idx_type b_nz = b.nnz (); 1359 octave_idx_type b_nz = b.nnz ();
1360 retval = SparseComplexMatrix (nc, b_nc, b_nz); 1360 retval = SparseComplexMatrix (nc, b_nc, b_nz);
1361 1361
1362 retval.xcidx(0) = 0; 1362 retval.xcidx(0) = 0;
1363 octave_idx_type ii = 0; 1363 octave_idx_type ii = 0;
1364 if (typ == SparseType::Diagonal) 1364 if (typ == MatrixType::Diagonal)
1365 for (octave_idx_type j = 0; j < b.cols(); j++) 1365 for (octave_idx_type j = 0; j < b.cols(); j++)
1366 { 1366 {
1367 for (octave_idx_type i = b.cidx(j); i < b.cidx(j+1); i++) 1367 for (octave_idx_type i = b.cidx(j); i < b.cidx(j+1); i++)
1368 { 1368 {
1369 if (b.ridx(i) >= nm) 1369 if (b.ridx(i) >= nm)
1418 1418
1419 return retval; 1419 return retval;
1420 } 1420 }
1421 1421
1422 ComplexMatrix 1422 ComplexMatrix
1423 SparseComplexMatrix::utsolve (SparseType &mattype, const Matrix& b, 1423 SparseComplexMatrix::utsolve (MatrixType &mattype, const Matrix& b,
1424 octave_idx_type& err, double& rcond, 1424 octave_idx_type& err, double& rcond,
1425 solve_singularity_handler sing_handler, 1425 solve_singularity_handler sing_handler,
1426 bool calc_cond) const 1426 bool calc_cond) const
1427 { 1427 {
1428 ComplexMatrix retval; 1428 ComplexMatrix retval;
1439 { 1439 {
1440 // Print spparms("spumoni") info if requested 1440 // Print spparms("spumoni") info if requested
1441 int typ = mattype.type (); 1441 int typ = mattype.type ();
1442 mattype.info (); 1442 mattype.info ();
1443 1443
1444 if (typ == SparseType::Permuted_Upper || 1444 if (typ == MatrixType::Permuted_Upper ||
1445 typ == SparseType::Upper) 1445 typ == MatrixType::Upper)
1446 { 1446 {
1447 double anorm = 0.; 1447 double anorm = 0.;
1448 double ainvnorm = 0.; 1448 double ainvnorm = 0.;
1449 octave_idx_type b_nc = b.cols (); 1449 octave_idx_type b_nc = b.cols ();
1450 rcond = 1.; 1450 rcond = 1.;
1460 if (atmp > anorm) 1460 if (atmp > anorm)
1461 anorm = atmp; 1461 anorm = atmp;
1462 } 1462 }
1463 } 1463 }
1464 1464
1465 if (typ == SparseType::Permuted_Upper) 1465 if (typ == MatrixType::Permuted_Upper)
1466 { 1466 {
1467 retval.resize (nc, b_nc); 1467 retval.resize (nc, b_nc);
1468 octave_idx_type *perm = mattype.triangular_perm (); 1468 octave_idx_type *perm = mattype.triangular_perm ();
1469 OCTAVE_LOCAL_BUFFER (Complex, work, nm); 1469 OCTAVE_LOCAL_BUFFER (Complex, work, nm);
1470 1470
1652 1652
1653 return retval; 1653 return retval;
1654 } 1654 }
1655 1655
1656 SparseComplexMatrix 1656 SparseComplexMatrix
1657 SparseComplexMatrix::utsolve (SparseType &mattype, const SparseMatrix& b, 1657 SparseComplexMatrix::utsolve (MatrixType &mattype, const SparseMatrix& b,
1658 octave_idx_type& err, double& rcond, 1658 octave_idx_type& err, double& rcond,
1659 solve_singularity_handler sing_handler, 1659 solve_singularity_handler sing_handler,
1660 bool calc_cond) const 1660 bool calc_cond) const
1661 { 1661 {
1662 SparseComplexMatrix retval; 1662 SparseComplexMatrix retval;
1673 { 1673 {
1674 // Print spparms("spumoni") info if requested 1674 // Print spparms("spumoni") info if requested
1675 int typ = mattype.type (); 1675 int typ = mattype.type ();
1676 mattype.info (); 1676 mattype.info ();
1677 1677
1678 if (typ == SparseType::Permuted_Upper || 1678 if (typ == MatrixType::Permuted_Upper ||
1679 typ == SparseType::Upper) 1679 typ == MatrixType::Upper)
1680 { 1680 {
1681 double anorm = 0.; 1681 double anorm = 0.;
1682 double ainvnorm = 0.; 1682 double ainvnorm = 0.;
1683 rcond = 1.; 1683 rcond = 1.;
1684 1684
1700 retval = SparseComplexMatrix (nc, b_nc, b_nz); 1700 retval = SparseComplexMatrix (nc, b_nc, b_nz);
1701 retval.xcidx(0) = 0; 1701 retval.xcidx(0) = 0;
1702 octave_idx_type ii = 0; 1702 octave_idx_type ii = 0;
1703 octave_idx_type x_nz = b_nz; 1703 octave_idx_type x_nz = b_nz;
1704 1704
1705 if (typ == SparseType::Permuted_Upper) 1705 if (typ == MatrixType::Permuted_Upper)
1706 { 1706 {
1707 octave_idx_type *perm = mattype.triangular_perm (); 1707 octave_idx_type *perm = mattype.triangular_perm ();
1708 OCTAVE_LOCAL_BUFFER (Complex, work, nm); 1708 OCTAVE_LOCAL_BUFFER (Complex, work, nm);
1709 1709
1710 OCTAVE_LOCAL_BUFFER (octave_idx_type, rperm, nc); 1710 OCTAVE_LOCAL_BUFFER (octave_idx_type, rperm, nc);
1937 } 1937 }
1938 return retval; 1938 return retval;
1939 } 1939 }
1940 1940
1941 ComplexMatrix 1941 ComplexMatrix
1942 SparseComplexMatrix::utsolve (SparseType &mattype, const ComplexMatrix& b, 1942 SparseComplexMatrix::utsolve (MatrixType &mattype, const ComplexMatrix& b,
1943 octave_idx_type& err, double& rcond, 1943 octave_idx_type& err, double& rcond,
1944 solve_singularity_handler sing_handler, 1944 solve_singularity_handler sing_handler,
1945 bool calc_cond) const 1945 bool calc_cond) const
1946 { 1946 {
1947 ComplexMatrix retval; 1947 ComplexMatrix retval;
1958 { 1958 {
1959 // Print spparms("spumoni") info if requested 1959 // Print spparms("spumoni") info if requested
1960 int typ = mattype.type (); 1960 int typ = mattype.type ();
1961 mattype.info (); 1961 mattype.info ();
1962 1962
1963 if (typ == SparseType::Permuted_Upper || 1963 if (typ == MatrixType::Permuted_Upper ||
1964 typ == SparseType::Upper) 1964 typ == MatrixType::Upper)
1965 { 1965 {
1966 double anorm = 0.; 1966 double anorm = 0.;
1967 double ainvnorm = 0.; 1967 double ainvnorm = 0.;
1968 octave_idx_type b_nc = b.cols (); 1968 octave_idx_type b_nc = b.cols ();
1969 rcond = 1.; 1969 rcond = 1.;
1979 if (atmp > anorm) 1979 if (atmp > anorm)
1980 anorm = atmp; 1980 anorm = atmp;
1981 } 1981 }
1982 } 1982 }
1983 1983
1984 if (typ == SparseType::Permuted_Upper) 1984 if (typ == MatrixType::Permuted_Upper)
1985 { 1985 {
1986 retval.resize (nc, b_nc); 1986 retval.resize (nc, b_nc);
1987 octave_idx_type *perm = mattype.triangular_perm (); 1987 octave_idx_type *perm = mattype.triangular_perm ();
1988 OCTAVE_LOCAL_BUFFER (Complex, work, nm); 1988 OCTAVE_LOCAL_BUFFER (Complex, work, nm);
1989 1989
2171 2171
2172 return retval; 2172 return retval;
2173 } 2173 }
2174 2174
2175 SparseComplexMatrix 2175 SparseComplexMatrix
2176 SparseComplexMatrix::utsolve (SparseType &mattype, const SparseComplexMatrix& b, 2176 SparseComplexMatrix::utsolve (MatrixType &mattype, const SparseComplexMatrix& b,
2177 octave_idx_type& err, double& rcond, 2177 octave_idx_type& err, double& rcond,
2178 solve_singularity_handler sing_handler, 2178 solve_singularity_handler sing_handler,
2179 bool calc_cond) const 2179 bool calc_cond) const
2180 { 2180 {
2181 SparseComplexMatrix retval; 2181 SparseComplexMatrix retval;
2192 { 2192 {
2193 // Print spparms("spumoni") info if requested 2193 // Print spparms("spumoni") info if requested
2194 int typ = mattype.type (); 2194 int typ = mattype.type ();
2195 mattype.info (); 2195 mattype.info ();
2196 2196
2197 if (typ == SparseType::Permuted_Upper || 2197 if (typ == MatrixType::Permuted_Upper ||
2198 typ == SparseType::Upper) 2198 typ == MatrixType::Upper)
2199 { 2199 {
2200 double anorm = 0.; 2200 double anorm = 0.;
2201 double ainvnorm = 0.; 2201 double ainvnorm = 0.;
2202 rcond = 1.; 2202 rcond = 1.;
2203 2203
2219 retval = SparseComplexMatrix (nc, b_nc, b_nz); 2219 retval = SparseComplexMatrix (nc, b_nc, b_nz);
2220 retval.xcidx(0) = 0; 2220 retval.xcidx(0) = 0;
2221 octave_idx_type ii = 0; 2221 octave_idx_type ii = 0;
2222 octave_idx_type x_nz = b_nz; 2222 octave_idx_type x_nz = b_nz;
2223 2223
2224 if (typ == SparseType::Permuted_Upper) 2224 if (typ == MatrixType::Permuted_Upper)
2225 { 2225 {
2226 octave_idx_type *perm = mattype.triangular_perm (); 2226 octave_idx_type *perm = mattype.triangular_perm ();
2227 OCTAVE_LOCAL_BUFFER (Complex, work, nm); 2227 OCTAVE_LOCAL_BUFFER (Complex, work, nm);
2228 2228
2229 OCTAVE_LOCAL_BUFFER (octave_idx_type, rperm, nc); 2229 OCTAVE_LOCAL_BUFFER (octave_idx_type, rperm, nc);
2457 2457
2458 return retval; 2458 return retval;
2459 } 2459 }
2460 2460
2461 ComplexMatrix 2461 ComplexMatrix
2462 SparseComplexMatrix::ltsolve (SparseType &mattype, const Matrix& b, 2462 SparseComplexMatrix::ltsolve (MatrixType &mattype, const Matrix& b,
2463 octave_idx_type& err, double& rcond, 2463 octave_idx_type& err, double& rcond,
2464 solve_singularity_handler sing_handler, 2464 solve_singularity_handler sing_handler,
2465 bool calc_cond) const 2465 bool calc_cond) const
2466 { 2466 {
2467 ComplexMatrix retval; 2467 ComplexMatrix retval;
2478 { 2478 {
2479 // Print spparms("spumoni") info if requested 2479 // Print spparms("spumoni") info if requested
2480 int typ = mattype.type (); 2480 int typ = mattype.type ();
2481 mattype.info (); 2481 mattype.info ();
2482 2482
2483 if (typ == SparseType::Permuted_Lower || 2483 if (typ == MatrixType::Permuted_Lower ||
2484 typ == SparseType::Lower) 2484 typ == MatrixType::Lower)
2485 { 2485 {
2486 double anorm = 0.; 2486 double anorm = 0.;
2487 double ainvnorm = 0.; 2487 double ainvnorm = 0.;
2488 octave_idx_type b_nc = b.cols (); 2488 octave_idx_type b_nc = b.cols ();
2489 rcond = 1.; 2489 rcond = 1.;
2499 if (atmp > anorm) 2499 if (atmp > anorm)
2500 anorm = atmp; 2500 anorm = atmp;
2501 } 2501 }
2502 } 2502 }
2503 2503
2504 if (typ == SparseType::Permuted_Lower) 2504 if (typ == MatrixType::Permuted_Lower)
2505 { 2505 {
2506 retval.resize (nc, b_nc); 2506 retval.resize (nc, b_nc);
2507 OCTAVE_LOCAL_BUFFER (Complex, work, nm); 2507 OCTAVE_LOCAL_BUFFER (Complex, work, nm);
2508 octave_idx_type *perm = mattype.triangular_perm (); 2508 octave_idx_type *perm = mattype.triangular_perm ();
2509 2509
2711 2711
2712 return retval; 2712 return retval;
2713 } 2713 }
2714 2714
2715 SparseComplexMatrix 2715 SparseComplexMatrix
2716 SparseComplexMatrix::ltsolve (SparseType &mattype, const SparseMatrix& b, 2716 SparseComplexMatrix::ltsolve (MatrixType &mattype, const SparseMatrix& b,
2717 octave_idx_type& err, double& rcond, 2717 octave_idx_type& err, double& rcond,
2718 solve_singularity_handler sing_handler, 2718 solve_singularity_handler sing_handler,
2719 bool calc_cond) const 2719 bool calc_cond) const
2720 { 2720 {
2721 SparseComplexMatrix retval; 2721 SparseComplexMatrix retval;
2733 { 2733 {
2734 // Print spparms("spumoni") info if requested 2734 // Print spparms("spumoni") info if requested
2735 int typ = mattype.type (); 2735 int typ = mattype.type ();
2736 mattype.info (); 2736 mattype.info ();
2737 2737
2738 if (typ == SparseType::Permuted_Lower || 2738 if (typ == MatrixType::Permuted_Lower ||
2739 typ == SparseType::Lower) 2739 typ == MatrixType::Lower)
2740 { 2740 {
2741 double anorm = 0.; 2741 double anorm = 0.;
2742 double ainvnorm = 0.; 2742 double ainvnorm = 0.;
2743 rcond = 1.; 2743 rcond = 1.;
2744 2744
2760 retval = SparseComplexMatrix (nc, b_nc, b_nz); 2760 retval = SparseComplexMatrix (nc, b_nc, b_nz);
2761 retval.xcidx(0) = 0; 2761 retval.xcidx(0) = 0;
2762 octave_idx_type ii = 0; 2762 octave_idx_type ii = 0;
2763 octave_idx_type x_nz = b_nz; 2763 octave_idx_type x_nz = b_nz;
2764 2764
2765 if (typ == SparseType::Permuted_Lower) 2765 if (typ == MatrixType::Permuted_Lower)
2766 { 2766 {
2767 OCTAVE_LOCAL_BUFFER (Complex, work, nm); 2767 OCTAVE_LOCAL_BUFFER (Complex, work, nm);
2768 octave_idx_type *perm = mattype.triangular_perm (); 2768 octave_idx_type *perm = mattype.triangular_perm ();
2769 2769
2770 for (octave_idx_type j = 0; j < b_nc; j++) 2770 for (octave_idx_type j = 0; j < b_nc; j++)
3017 3017
3018 return retval; 3018 return retval;
3019 } 3019 }
3020 3020
3021 ComplexMatrix 3021 ComplexMatrix
3022 SparseComplexMatrix::ltsolve (SparseType &mattype, const ComplexMatrix& b, 3022 SparseComplexMatrix::ltsolve (MatrixType &mattype, const ComplexMatrix& b,
3023 octave_idx_type& err, double& rcond, 3023 octave_idx_type& err, double& rcond,
3024 solve_singularity_handler sing_handler, 3024 solve_singularity_handler sing_handler,
3025 bool calc_cond) const 3025 bool calc_cond) const
3026 { 3026 {
3027 ComplexMatrix retval; 3027 ComplexMatrix retval;
3038 { 3038 {
3039 // Print spparms("spumoni") info if requested 3039 // Print spparms("spumoni") info if requested
3040 int typ = mattype.type (); 3040 int typ = mattype.type ();
3041 mattype.info (); 3041 mattype.info ();
3042 3042
3043 if (typ == SparseType::Permuted_Lower || 3043 if (typ == MatrixType::Permuted_Lower ||
3044 typ == SparseType::Lower) 3044 typ == MatrixType::Lower)
3045 { 3045 {
3046 double anorm = 0.; 3046 double anorm = 0.;
3047 double ainvnorm = 0.; 3047 double ainvnorm = 0.;
3048 octave_idx_type b_nc = b.cols (); 3048 octave_idx_type b_nc = b.cols ();
3049 rcond = 1.; 3049 rcond = 1.;
3059 if (atmp > anorm) 3059 if (atmp > anorm)
3060 anorm = atmp; 3060 anorm = atmp;
3061 } 3061 }
3062 } 3062 }
3063 3063
3064 if (typ == SparseType::Permuted_Lower) 3064 if (typ == MatrixType::Permuted_Lower)
3065 { 3065 {
3066 retval.resize (nc, b_nc); 3066 retval.resize (nc, b_nc);
3067 OCTAVE_LOCAL_BUFFER (Complex, work, nm); 3067 OCTAVE_LOCAL_BUFFER (Complex, work, nm);
3068 octave_idx_type *perm = mattype.triangular_perm (); 3068 octave_idx_type *perm = mattype.triangular_perm ();
3069 3069
3275 3275
3276 return retval; 3276 return retval;
3277 } 3277 }
3278 3278
3279 SparseComplexMatrix 3279 SparseComplexMatrix
3280 SparseComplexMatrix::ltsolve (SparseType &mattype, const SparseComplexMatrix& b, 3280 SparseComplexMatrix::ltsolve (MatrixType &mattype, const SparseComplexMatrix& b,
3281 octave_idx_type& err, double& rcond, 3281 octave_idx_type& err, double& rcond,
3282 solve_singularity_handler sing_handler, 3282 solve_singularity_handler sing_handler,
3283 bool calc_cond) const 3283 bool calc_cond) const
3284 { 3284 {
3285 SparseComplexMatrix retval; 3285 SparseComplexMatrix retval;
3296 { 3296 {
3297 // Print spparms("spumoni") info if requested 3297 // Print spparms("spumoni") info if requested
3298 int typ = mattype.type (); 3298 int typ = mattype.type ();
3299 mattype.info (); 3299 mattype.info ();
3300 3300
3301 if (typ == SparseType::Permuted_Lower || 3301 if (typ == MatrixType::Permuted_Lower ||
3302 typ == SparseType::Lower) 3302 typ == MatrixType::Lower)
3303 { 3303 {
3304 double anorm = 0.; 3304 double anorm = 0.;
3305 double ainvnorm = 0.; 3305 double ainvnorm = 0.;
3306 rcond = 1.; 3306 rcond = 1.;
3307 3307
3323 retval = SparseComplexMatrix (nc, b_nc, b_nz); 3323 retval = SparseComplexMatrix (nc, b_nc, b_nz);
3324 retval.xcidx(0) = 0; 3324 retval.xcidx(0) = 0;
3325 octave_idx_type ii = 0; 3325 octave_idx_type ii = 0;
3326 octave_idx_type x_nz = b_nz; 3326 octave_idx_type x_nz = b_nz;
3327 3327
3328 if (typ == SparseType::Permuted_Lower) 3328 if (typ == MatrixType::Permuted_Lower)
3329 { 3329 {
3330 OCTAVE_LOCAL_BUFFER (Complex, work, nm); 3330 OCTAVE_LOCAL_BUFFER (Complex, work, nm);
3331 octave_idx_type *perm = mattype.triangular_perm (); 3331 octave_idx_type *perm = mattype.triangular_perm ();
3332 3332
3333 for (octave_idx_type j = 0; j < b_nc; j++) 3333 for (octave_idx_type j = 0; j < b_nc; j++)
3580 3580
3581 return retval; 3581 return retval;
3582 } 3582 }
3583 3583
3584 ComplexMatrix 3584 ComplexMatrix
3585 SparseComplexMatrix::trisolve (SparseType &mattype, const Matrix& b, 3585 SparseComplexMatrix::trisolve (MatrixType &mattype, const Matrix& b,
3586 octave_idx_type& err, double& rcond, 3586 octave_idx_type& err, double& rcond,
3587 solve_singularity_handler sing_handler, 3587 solve_singularity_handler sing_handler,
3588 bool calc_cond) const 3588 bool calc_cond) const
3589 { 3589 {
3590 ComplexMatrix retval; 3590 ComplexMatrix retval;
3603 { 3603 {
3604 // Print spparms("spumoni") info if requested 3604 // Print spparms("spumoni") info if requested
3605 volatile int typ = mattype.type (); 3605 volatile int typ = mattype.type ();
3606 mattype.info (); 3606 mattype.info ();
3607 3607
3608 if (typ == SparseType::Tridiagonal_Hermitian) 3608 if (typ == MatrixType::Tridiagonal_Hermitian)
3609 { 3609 {
3610 OCTAVE_LOCAL_BUFFER (double, D, nr); 3610 OCTAVE_LOCAL_BUFFER (double, D, nr);
3611 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1); 3611 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1);
3612 3612
3613 if (mattype.is_dense ()) 3613 if (mattype.is_dense ())
3653 ("unrecoverable error in zptsv"); 3653 ("unrecoverable error in zptsv");
3654 else if (err != 0) 3654 else if (err != 0)
3655 { 3655 {
3656 err = 0; 3656 err = 0;
3657 mattype.mark_as_unsymmetric (); 3657 mattype.mark_as_unsymmetric ();
3658 typ = SparseType::Tridiagonal; 3658 typ = MatrixType::Tridiagonal;
3659 } 3659 }
3660 else 3660 else
3661 rcond = 1.; 3661 rcond = 1.;
3662 } 3662 }
3663 3663
3664 if (typ == SparseType::Tridiagonal) 3664 if (typ == MatrixType::Tridiagonal)
3665 { 3665 {
3666 OCTAVE_LOCAL_BUFFER (Complex, DU, nr - 1); 3666 OCTAVE_LOCAL_BUFFER (Complex, DU, nr - 1);
3667 OCTAVE_LOCAL_BUFFER (Complex, D, nr); 3667 OCTAVE_LOCAL_BUFFER (Complex, D, nr);
3668 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1); 3668 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1);
3669 3669
3727 3727
3728 } 3728 }
3729 else 3729 else
3730 rcond = 1.; 3730 rcond = 1.;
3731 } 3731 }
3732 else if (typ != SparseType::Tridiagonal_Hermitian) 3732 else if (typ != MatrixType::Tridiagonal_Hermitian)
3733 (*current_liboctave_error_handler) ("incorrect matrix type"); 3733 (*current_liboctave_error_handler) ("incorrect matrix type");
3734 } 3734 }
3735 3735
3736 return retval; 3736 return retval;
3737 } 3737 }
3738 3738
3739 SparseComplexMatrix 3739 SparseComplexMatrix
3740 SparseComplexMatrix::trisolve (SparseType &mattype, const SparseMatrix& b, 3740 SparseComplexMatrix::trisolve (MatrixType &mattype, const SparseMatrix& b,
3741 octave_idx_type& err, double& rcond, 3741 octave_idx_type& err, double& rcond,
3742 solve_singularity_handler sing_handler, 3742 solve_singularity_handler sing_handler,
3743 bool calc_cond) const 3743 bool calc_cond) const
3744 { 3744 {
3745 SparseComplexMatrix retval; 3745 SparseComplexMatrix retval;
3759 // Print spparms("spumoni") info if requested 3759 // Print spparms("spumoni") info if requested
3760 int typ = mattype.type (); 3760 int typ = mattype.type ();
3761 mattype.info (); 3761 mattype.info ();
3762 3762
3763 // Note can't treat symmetric case as there is no dpttrf function 3763 // Note can't treat symmetric case as there is no dpttrf function
3764 if (typ == SparseType::Tridiagonal || 3764 if (typ == MatrixType::Tridiagonal ||
3765 typ == SparseType::Tridiagonal_Hermitian) 3765 typ == MatrixType::Tridiagonal_Hermitian)
3766 { 3766 {
3767 OCTAVE_LOCAL_BUFFER (Complex, DU2, nr - 2); 3767 OCTAVE_LOCAL_BUFFER (Complex, DU2, nr - 2);
3768 OCTAVE_LOCAL_BUFFER (Complex, DU, nr - 1); 3768 OCTAVE_LOCAL_BUFFER (Complex, DU, nr - 1);
3769 OCTAVE_LOCAL_BUFFER (Complex, D, nr); 3769 OCTAVE_LOCAL_BUFFER (Complex, D, nr);
3770 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1); 3770 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1);
3885 3885
3886 retval.maybe_compress (); 3886 retval.maybe_compress ();
3887 } 3887 }
3888 } 3888 }
3889 } 3889 }
3890 else if (typ != SparseType::Tridiagonal_Hermitian) 3890 else if (typ != MatrixType::Tridiagonal_Hermitian)
3891 (*current_liboctave_error_handler) ("incorrect matrix type"); 3891 (*current_liboctave_error_handler) ("incorrect matrix type");
3892 } 3892 }
3893 3893
3894 return retval; 3894 return retval;
3895 } 3895 }
3896 3896
3897 ComplexMatrix 3897 ComplexMatrix
3898 SparseComplexMatrix::trisolve (SparseType &mattype, const ComplexMatrix& b, 3898 SparseComplexMatrix::trisolve (MatrixType &mattype, const ComplexMatrix& b,
3899 octave_idx_type& err, double& rcond, 3899 octave_idx_type& err, double& rcond,
3900 solve_singularity_handler sing_handler, 3900 solve_singularity_handler sing_handler,
3901 bool calc_cond) const 3901 bool calc_cond) const
3902 { 3902 {
3903 ComplexMatrix retval; 3903 ComplexMatrix retval;
3916 { 3916 {
3917 // Print spparms("spumoni") info if requested 3917 // Print spparms("spumoni") info if requested
3918 volatile int typ = mattype.type (); 3918 volatile int typ = mattype.type ();
3919 mattype.info (); 3919 mattype.info ();
3920 3920
3921 if (typ == SparseType::Tridiagonal_Hermitian) 3921 if (typ == MatrixType::Tridiagonal_Hermitian)
3922 { 3922 {
3923 OCTAVE_LOCAL_BUFFER (double, D, nr); 3923 OCTAVE_LOCAL_BUFFER (double, D, nr);
3924 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1); 3924 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1);
3925 3925
3926 if (mattype.is_dense ()) 3926 if (mattype.is_dense ())
3972 } 3972 }
3973 else if (err != 0) 3973 else if (err != 0)
3974 { 3974 {
3975 err = 0; 3975 err = 0;
3976 mattype.mark_as_unsymmetric (); 3976 mattype.mark_as_unsymmetric ();
3977 typ = SparseType::Tridiagonal; 3977 typ = MatrixType::Tridiagonal;
3978 } 3978 }
3979 } 3979 }
3980 3980
3981 if (typ == SparseType::Tridiagonal) 3981 if (typ == MatrixType::Tridiagonal)
3982 { 3982 {
3983 OCTAVE_LOCAL_BUFFER (Complex, DU, nr - 1); 3983 OCTAVE_LOCAL_BUFFER (Complex, DU, nr - 1);
3984 OCTAVE_LOCAL_BUFFER (Complex, D, nr); 3984 OCTAVE_LOCAL_BUFFER (Complex, D, nr);
3985 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1); 3985 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1);
3986 3986
4047 else 4047 else
4048 (*current_liboctave_error_handler) 4048 (*current_liboctave_error_handler)
4049 ("matrix singular to machine precision"); 4049 ("matrix singular to machine precision");
4050 } 4050 }
4051 } 4051 }
4052 else if (typ != SparseType::Tridiagonal_Hermitian) 4052 else if (typ != MatrixType::Tridiagonal_Hermitian)
4053 (*current_liboctave_error_handler) ("incorrect matrix type"); 4053 (*current_liboctave_error_handler) ("incorrect matrix type");
4054 } 4054 }
4055 4055
4056 return retval; 4056 return retval;
4057 } 4057 }
4058 4058
4059 SparseComplexMatrix 4059 SparseComplexMatrix
4060 SparseComplexMatrix::trisolve (SparseType &mattype, 4060 SparseComplexMatrix::trisolve (MatrixType &mattype,
4061 const SparseComplexMatrix& b, 4061 const SparseComplexMatrix& b,
4062 octave_idx_type& err, double& rcond, 4062 octave_idx_type& err, double& rcond,
4063 solve_singularity_handler sing_handler, 4063 solve_singularity_handler sing_handler,
4064 bool calc_cond) const 4064 bool calc_cond) const
4065 { 4065 {
4080 // Print spparms("spumoni") info if requested 4080 // Print spparms("spumoni") info if requested
4081 int typ = mattype.type (); 4081 int typ = mattype.type ();
4082 mattype.info (); 4082 mattype.info ();
4083 4083
4084 // Note can't treat symmetric case as there is no dpttrf function 4084 // Note can't treat symmetric case as there is no dpttrf function
4085 if (typ == SparseType::Tridiagonal || 4085 if (typ == MatrixType::Tridiagonal ||
4086 typ == SparseType::Tridiagonal_Hermitian) 4086 typ == MatrixType::Tridiagonal_Hermitian)
4087 { 4087 {
4088 OCTAVE_LOCAL_BUFFER (Complex, DU2, nr - 2); 4088 OCTAVE_LOCAL_BUFFER (Complex, DU2, nr - 2);
4089 OCTAVE_LOCAL_BUFFER (Complex, DU, nr - 1); 4089 OCTAVE_LOCAL_BUFFER (Complex, DU, nr - 1);
4090 OCTAVE_LOCAL_BUFFER (Complex, D, nr); 4090 OCTAVE_LOCAL_BUFFER (Complex, D, nr);
4091 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1); 4091 OCTAVE_LOCAL_BUFFER (Complex, DL, nr - 1);
4217 4217
4218 retval.maybe_compress (); 4218 retval.maybe_compress ();
4219 } 4219 }
4220 } 4220 }
4221 } 4221 }
4222 else if (typ != SparseType::Tridiagonal_Hermitian) 4222 else if (typ != MatrixType::Tridiagonal_Hermitian)
4223 (*current_liboctave_error_handler) ("incorrect matrix type"); 4223 (*current_liboctave_error_handler) ("incorrect matrix type");
4224 } 4224 }
4225 4225
4226 return retval; 4226 return retval;
4227 } 4227 }
4228 4228
4229 ComplexMatrix 4229 ComplexMatrix
4230 SparseComplexMatrix::bsolve (SparseType &mattype, const Matrix& b, 4230 SparseComplexMatrix::bsolve (MatrixType &mattype, const Matrix& b,
4231 octave_idx_type& err, double& rcond, 4231 octave_idx_type& err, double& rcond,
4232 solve_singularity_handler sing_handler, 4232 solve_singularity_handler sing_handler,
4233 bool calc_cond) const 4233 bool calc_cond) const
4234 { 4234 {
4235 ComplexMatrix retval; 4235 ComplexMatrix retval;
4245 { 4245 {
4246 // Print spparms("spumoni") info if requested 4246 // Print spparms("spumoni") info if requested
4247 volatile int typ = mattype.type (); 4247 volatile int typ = mattype.type ();
4248 mattype.info (); 4248 mattype.info ();
4249 4249
4250 if (typ == SparseType::Banded_Hermitian) 4250 if (typ == MatrixType::Banded_Hermitian)
4251 { 4251 {
4252 octave_idx_type n_lower = mattype.nlower (); 4252 octave_idx_type n_lower = mattype.nlower ();
4253 octave_idx_type ldm = n_lower + 1; 4253 octave_idx_type ldm = n_lower + 1;
4254 ComplexMatrix m_band (ldm, nc); 4254 ComplexMatrix m_band (ldm, nc);
4255 Complex *tmp_data = m_band.fortran_vec (); 4255 Complex *tmp_data = m_band.fortran_vec ();
4290 { 4290 {
4291 rcond = 0.0; 4291 rcond = 0.0;
4292 // Matrix is not positive definite!! Fall through to 4292 // Matrix is not positive definite!! Fall through to
4293 // unsymmetric banded solver. 4293 // unsymmetric banded solver.
4294 mattype.mark_as_unsymmetric (); 4294 mattype.mark_as_unsymmetric ();
4295 typ = SparseType::Banded; 4295 typ = MatrixType::Banded;
4296 err = 0; 4296 err = 0;
4297 } 4297 }
4298 else 4298 else
4299 { 4299 {
4300 if (calc_cond) 4300 if (calc_cond)
4363 } 4363 }
4364 } 4364 }
4365 } 4365 }
4366 } 4366 }
4367 4367
4368 if (typ == SparseType::Banded) 4368 if (typ == MatrixType::Banded)
4369 { 4369 {
4370 // Create the storage for the banded form of the sparse matrix 4370 // Create the storage for the banded form of the sparse matrix
4371 octave_idx_type n_upper = mattype.nupper (); 4371 octave_idx_type n_upper = mattype.nupper ();
4372 octave_idx_type n_lower = mattype.nlower (); 4372 octave_idx_type n_lower = mattype.nlower ();
4373 octave_idx_type ldm = n_upper + 2 * n_lower + 1; 4373 octave_idx_type ldm = n_upper + 2 * n_lower + 1;
4491 ("unrecoverable error in zgbtrs"); 4491 ("unrecoverable error in zgbtrs");
4492 } 4492 }
4493 } 4493 }
4494 } 4494 }
4495 } 4495 }
4496 else if (typ != SparseType::Banded_Hermitian) 4496 else if (typ != MatrixType::Banded_Hermitian)
4497 (*current_liboctave_error_handler) ("incorrect matrix type"); 4497 (*current_liboctave_error_handler) ("incorrect matrix type");
4498 } 4498 }
4499 4499
4500 return retval; 4500 return retval;
4501 } 4501 }
4502 4502
4503 SparseComplexMatrix 4503 SparseComplexMatrix
4504 SparseComplexMatrix::bsolve (SparseType &mattype, const SparseMatrix& b, 4504 SparseComplexMatrix::bsolve (MatrixType &mattype, const SparseMatrix& b,
4505 octave_idx_type& err, double& rcond, 4505 octave_idx_type& err, double& rcond,
4506 solve_singularity_handler sing_handler, 4506 solve_singularity_handler sing_handler,
4507 bool calc_cond) const 4507 bool calc_cond) const
4508 { 4508 {
4509 SparseComplexMatrix retval; 4509 SparseComplexMatrix retval;
4519 { 4519 {
4520 // Print spparms("spumoni") info if requested 4520 // Print spparms("spumoni") info if requested
4521 volatile int typ = mattype.type (); 4521 volatile int typ = mattype.type ();
4522 mattype.info (); 4522 mattype.info ();
4523 4523
4524 if (typ == SparseType::Banded_Hermitian) 4524 if (typ == MatrixType::Banded_Hermitian)
4525 { 4525 {
4526 octave_idx_type n_lower = mattype.nlower (); 4526 octave_idx_type n_lower = mattype.nlower ();
4527 octave_idx_type ldm = n_lower + 1; 4527 octave_idx_type ldm = n_lower + 1;
4528 4528
4529 ComplexMatrix m_band (ldm, nc); 4529 ComplexMatrix m_band (ldm, nc);
4563 { 4563 {
4564 if (err != 0) 4564 if (err != 0)
4565 { 4565 {
4566 rcond = 0.0; 4566 rcond = 0.0;
4567 mattype.mark_as_unsymmetric (); 4567 mattype.mark_as_unsymmetric ();
4568 typ = SparseType::Banded; 4568 typ = MatrixType::Banded;
4569 err = 0; 4569 err = 0;
4570 } 4570 }
4571 else 4571 else
4572 { 4572 {
4573 if (calc_cond) 4573 if (calc_cond)
4675 } 4675 }
4676 } 4676 }
4677 } 4677 }
4678 } 4678 }
4679 4679
4680 if (typ == SparseType::Banded) 4680 if (typ == MatrixType::Banded)
4681 { 4681 {
4682 // Create the storage for the banded form of the sparse matrix 4682 // Create the storage for the banded form of the sparse matrix
4683 octave_idx_type n_upper = mattype.nupper (); 4683 octave_idx_type n_upper = mattype.nupper ();
4684 octave_idx_type n_lower = mattype.nlower (); 4684 octave_idx_type n_lower = mattype.nlower ();
4685 octave_idx_type ldm = n_upper + 2 * n_lower + 1; 4685 octave_idx_type ldm = n_upper + 2 * n_lower + 1;
4842 retval.maybe_compress (); 4842 retval.maybe_compress ();
4843 } 4843 }
4844 } 4844 }
4845 } 4845 }
4846 } 4846 }
4847 else if (typ != SparseType::Banded_Hermitian) 4847 else if (typ != MatrixType::Banded_Hermitian)
4848 (*current_liboctave_error_handler) ("incorrect matrix type"); 4848 (*current_liboctave_error_handler) ("incorrect matrix type");
4849 } 4849 }
4850 4850
4851 return retval; 4851 return retval;
4852 } 4852 }
4853 4853
4854 ComplexMatrix 4854 ComplexMatrix
4855 SparseComplexMatrix::bsolve (SparseType &mattype, const ComplexMatrix& b, 4855 SparseComplexMatrix::bsolve (MatrixType &mattype, const ComplexMatrix& b,
4856 octave_idx_type& err, double& rcond, 4856 octave_idx_type& err, double& rcond,
4857 solve_singularity_handler sing_handler, 4857 solve_singularity_handler sing_handler,
4858 bool calc_cond) const 4858 bool calc_cond) const
4859 { 4859 {
4860 ComplexMatrix retval; 4860 ComplexMatrix retval;
4870 { 4870 {
4871 // Print spparms("spumoni") info if requested 4871 // Print spparms("spumoni") info if requested
4872 volatile int typ = mattype.type (); 4872 volatile int typ = mattype.type ();
4873 mattype.info (); 4873 mattype.info ();
4874 4874
4875 if (typ == SparseType::Banded_Hermitian) 4875 if (typ == MatrixType::Banded_Hermitian)
4876 { 4876 {
4877 octave_idx_type n_lower = mattype.nlower (); 4877 octave_idx_type n_lower = mattype.nlower ();
4878 octave_idx_type ldm = n_lower + 1; 4878 octave_idx_type ldm = n_lower + 1;
4879 4879
4880 ComplexMatrix m_band (ldm, nc); 4880 ComplexMatrix m_band (ldm, nc);
4916 { 4916 {
4917 // Matrix is not positive definite!! Fall through to 4917 // Matrix is not positive definite!! Fall through to
4918 // unsymmetric banded solver. 4918 // unsymmetric banded solver.
4919 rcond = 0.0; 4919 rcond = 0.0;
4920 mattype.mark_as_unsymmetric (); 4920 mattype.mark_as_unsymmetric ();
4921 typ = SparseType::Banded; 4921 typ = MatrixType::Banded;
4922 err = 0; 4922 err = 0;
4923 } 4923 }
4924 else 4924 else
4925 { 4925 {
4926 if (calc_cond) 4926 if (calc_cond)
4992 } 4992 }
4993 } 4993 }
4994 } 4994 }
4995 } 4995 }
4996 4996
4997 if (typ == SparseType::Banded) 4997 if (typ == MatrixType::Banded)
4998 { 4998 {
4999 // Create the storage for the banded form of the sparse matrix 4999 // Create the storage for the banded form of the sparse matrix
5000 octave_idx_type n_upper = mattype.nupper (); 5000 octave_idx_type n_upper = mattype.nupper ();
5001 octave_idx_type n_lower = mattype.nlower (); 5001 octave_idx_type n_lower = mattype.nlower ();
5002 octave_idx_type ldm = n_upper + 2 * n_lower + 1; 5002 octave_idx_type ldm = n_upper + 2 * n_lower + 1;
5119 } 5119 }
5120 } 5120 }
5121 } 5121 }
5122 } 5122 }
5123 } 5123 }
5124 else if (typ != SparseType::Banded_Hermitian) 5124 else if (typ != MatrixType::Banded_Hermitian)
5125 (*current_liboctave_error_handler) ("incorrect matrix type"); 5125 (*current_liboctave_error_handler) ("incorrect matrix type");
5126 } 5126 }
5127 5127
5128 return retval; 5128 return retval;
5129 } 5129 }
5130 5130
5131 SparseComplexMatrix 5131 SparseComplexMatrix
5132 SparseComplexMatrix::bsolve (SparseType &mattype, const SparseComplexMatrix& b, 5132 SparseComplexMatrix::bsolve (MatrixType &mattype, const SparseComplexMatrix& b,
5133 octave_idx_type& err, double& rcond, 5133 octave_idx_type& err, double& rcond,
5134 solve_singularity_handler sing_handler, 5134 solve_singularity_handler sing_handler,
5135 bool calc_cond) const 5135 bool calc_cond) const
5136 { 5136 {
5137 SparseComplexMatrix retval; 5137 SparseComplexMatrix retval;
5147 { 5147 {
5148 // Print spparms("spumoni") info if requested 5148 // Print spparms("spumoni") info if requested
5149 volatile int typ = mattype.type (); 5149 volatile int typ = mattype.type ();
5150 mattype.info (); 5150 mattype.info ();
5151 5151
5152 if (typ == SparseType::Banded_Hermitian) 5152 if (typ == MatrixType::Banded_Hermitian)
5153 { 5153 {
5154 octave_idx_type n_lower = mattype.nlower (); 5154 octave_idx_type n_lower = mattype.nlower ();
5155 octave_idx_type ldm = n_lower + 1; 5155 octave_idx_type ldm = n_lower + 1;
5156 5156
5157 ComplexMatrix m_band (ldm, nc); 5157 ComplexMatrix m_band (ldm, nc);
5192 if (err != 0) 5192 if (err != 0)
5193 { 5193 {
5194 // Matrix is not positive definite!! Fall through to 5194 // Matrix is not positive definite!! Fall through to
5195 // unsymmetric banded solver. 5195 // unsymmetric banded solver.
5196 mattype.mark_as_unsymmetric (); 5196 mattype.mark_as_unsymmetric ();
5197 typ = SparseType::Banded; 5197 typ = MatrixType::Banded;
5198 5198
5199 rcond = 0.0; 5199 rcond = 0.0;
5200 err = 0; 5200 err = 0;
5201 } 5201 }
5202 else 5202 else
5311 } 5311 }
5312 } 5312 }
5313 } 5313 }
5314 } 5314 }
5315 5315
5316 if (typ == SparseType::Banded) 5316 if (typ == MatrixType::Banded)
5317 { 5317 {
5318 // Create the storage for the banded form of the sparse matrix 5318 // Create the storage for the banded form of the sparse matrix
5319 octave_idx_type n_upper = mattype.nupper (); 5319 octave_idx_type n_upper = mattype.nupper ();
5320 octave_idx_type n_lower = mattype.nlower (); 5320 octave_idx_type n_lower = mattype.nlower ();
5321 octave_idx_type ldm = n_upper + 2 * n_lower + 1; 5321 octave_idx_type ldm = n_upper + 2 * n_lower + 1;
5479 retval.maybe_compress (); 5479 retval.maybe_compress ();
5480 } 5480 }
5481 } 5481 }
5482 } 5482 }
5483 } 5483 }
5484 else if (typ != SparseType::Banded_Hermitian) 5484 else if (typ != MatrixType::Banded_Hermitian)
5485 (*current_liboctave_error_handler) ("incorrect matrix type"); 5485 (*current_liboctave_error_handler) ("incorrect matrix type");
5486 } 5486 }
5487 5487
5488 return retval; 5488 return retval;
5489 } 5489 }
5603 5603
5604 return Numeric; 5604 return Numeric;
5605 } 5605 }
5606 5606
5607 ComplexMatrix 5607 ComplexMatrix
5608 SparseComplexMatrix::fsolve (SparseType &mattype, const Matrix& b, 5608 SparseComplexMatrix::fsolve (MatrixType &mattype, const Matrix& b,
5609 octave_idx_type& err, double& rcond, 5609 octave_idx_type& err, double& rcond,
5610 solve_singularity_handler sing_handler, 5610 solve_singularity_handler sing_handler,
5611 bool calc_cond) const 5611 bool calc_cond) const
5612 { 5612 {
5613 ComplexMatrix retval; 5613 ComplexMatrix retval;
5623 { 5623 {
5624 // Print spparms("spumoni") info if requested 5624 // Print spparms("spumoni") info if requested
5625 volatile int typ = mattype.type (); 5625 volatile int typ = mattype.type ();
5626 mattype.info (); 5626 mattype.info ();
5627 5627
5628 if (typ == SparseType::Hermitian) 5628 if (typ == MatrixType::Hermitian)
5629 { 5629 {
5630 #ifdef HAVE_CHOLMOD 5630 #ifdef HAVE_CHOLMOD
5631 cholmod_common Common; 5631 cholmod_common Common;
5632 cholmod_common *cm = &Common; 5632 cholmod_common *cm = &Common;
5633 5633
5722 5722
5723 if (rcond == 0.0) 5723 if (rcond == 0.0)
5724 { 5724 {
5725 // Either its indefinite or singular. Try UMFPACK 5725 // Either its indefinite or singular. Try UMFPACK
5726 mattype.mark_as_unsymmetric (); 5726 mattype.mark_as_unsymmetric ();
5727 typ = SparseType::Full; 5727 typ = MatrixType::Full;
5728 } 5728 }
5729 else 5729 else
5730 { 5730 {
5731 volatile double rcond_plus_one = rcond + 1.0; 5731 volatile double rcond_plus_one = rcond + 1.0;
5732 5732
5770 #else 5770 #else
5771 (*current_liboctave_warning_handler) 5771 (*current_liboctave_warning_handler)
5772 ("CHOLMOD not installed"); 5772 ("CHOLMOD not installed");
5773 5773
5774 mattype.mark_as_unsymmetric (); 5774 mattype.mark_as_unsymmetric ();
5775 typ = SparseType::Full; 5775 typ = MatrixType::Full;
5776 #endif 5776 #endif
5777 } 5777 }
5778 5778
5779 if (typ == SparseType::Full) 5779 if (typ == MatrixType::Full)
5780 { 5780 {
5781 #ifdef HAVE_UMFPACK 5781 #ifdef HAVE_UMFPACK
5782 Matrix Control, Info; 5782 Matrix Control, Info;
5783 void *Numeric = factorize (err, rcond, Control, Info, 5783 void *Numeric = factorize (err, rcond, Control, Info,
5784 sing_handler, calc_cond); 5784 sing_handler, calc_cond);
5852 5852
5853 #else 5853 #else
5854 (*current_liboctave_error_handler) ("UMFPACK not installed"); 5854 (*current_liboctave_error_handler) ("UMFPACK not installed");
5855 #endif 5855 #endif
5856 } 5856 }
5857 else if (typ != SparseType::Hermitian) 5857 else if (typ != MatrixType::Hermitian)
5858 (*current_liboctave_error_handler) ("incorrect matrix type"); 5858 (*current_liboctave_error_handler) ("incorrect matrix type");
5859 } 5859 }
5860 5860
5861 return retval; 5861 return retval;
5862 } 5862 }
5863 5863
5864 SparseComplexMatrix 5864 SparseComplexMatrix
5865 SparseComplexMatrix::fsolve (SparseType &mattype, const SparseMatrix& b, 5865 SparseComplexMatrix::fsolve (MatrixType &mattype, const SparseMatrix& b,
5866 octave_idx_type& err, double& rcond, 5866 octave_idx_type& err, double& rcond,
5867 solve_singularity_handler sing_handler, 5867 solve_singularity_handler sing_handler,
5868 bool calc_cond) const 5868 bool calc_cond) const
5869 { 5869 {
5870 SparseComplexMatrix retval; 5870 SparseComplexMatrix retval;
5880 { 5880 {
5881 // Print spparms("spumoni") info if requested 5881 // Print spparms("spumoni") info if requested
5882 volatile int typ = mattype.type (); 5882 volatile int typ = mattype.type ();
5883 mattype.info (); 5883 mattype.info ();
5884 5884
5885 if (typ == SparseType::Hermitian) 5885 if (typ == MatrixType::Hermitian)
5886 { 5886 {
5887 #ifdef HAVE_CHOLMOD 5887 #ifdef HAVE_CHOLMOD
5888 cholmod_common Common; 5888 cholmod_common Common;
5889 cholmod_common *cm = &Common; 5889 cholmod_common *cm = &Common;
5890 5890
5992 5992
5993 if (rcond == 0.0) 5993 if (rcond == 0.0)
5994 { 5994 {
5995 // Either its indefinite or singular. Try UMFPACK 5995 // Either its indefinite or singular. Try UMFPACK
5996 mattype.mark_as_unsymmetric (); 5996 mattype.mark_as_unsymmetric ();
5997 typ = SparseType::Full; 5997 typ = MatrixType::Full;
5998 } 5998 }
5999 else 5999 else
6000 { 6000 {
6001 volatile double rcond_plus_one = rcond + 1.0; 6001 volatile double rcond_plus_one = rcond + 1.0;
6002 6002
6046 #else 6046 #else
6047 (*current_liboctave_warning_handler) 6047 (*current_liboctave_warning_handler)
6048 ("CHOLMOD not installed"); 6048 ("CHOLMOD not installed");
6049 6049
6050 mattype.mark_as_unsymmetric (); 6050 mattype.mark_as_unsymmetric ();
6051 typ = SparseType::Full; 6051 typ = MatrixType::Full;
6052 #endif 6052 #endif
6053 } 6053 }
6054 6054
6055 if (typ == SparseType::Full) 6055 if (typ == MatrixType::Full)
6056 { 6056 {
6057 #ifdef HAVE_UMFPACK 6057 #ifdef HAVE_UMFPACK
6058 Matrix Control, Info; 6058 Matrix Control, Info;
6059 void *Numeric = factorize (err, rcond, Control, Info, 6059 void *Numeric = factorize (err, rcond, Control, Info,
6060 sing_handler, calc_cond); 6060 sing_handler, calc_cond);
6160 6160
6161 #else 6161 #else
6162 (*current_liboctave_error_handler) ("UMFPACK not installed"); 6162 (*current_liboctave_error_handler) ("UMFPACK not installed");
6163 #endif 6163 #endif
6164 } 6164 }
6165 else if (typ != SparseType::Hermitian) 6165 else if (typ != MatrixType::Hermitian)
6166 (*current_liboctave_error_handler) ("incorrect matrix type"); 6166 (*current_liboctave_error_handler) ("incorrect matrix type");
6167 } 6167 }
6168 6168
6169 return retval; 6169 return retval;
6170 } 6170 }
6171 6171
6172 ComplexMatrix 6172 ComplexMatrix
6173 SparseComplexMatrix::fsolve (SparseType &mattype, const ComplexMatrix& b, 6173 SparseComplexMatrix::fsolve (MatrixType &mattype, const ComplexMatrix& b,
6174 octave_idx_type& err, double& rcond, 6174 octave_idx_type& err, double& rcond,
6175 solve_singularity_handler sing_handler, 6175 solve_singularity_handler sing_handler,
6176 bool calc_cond) const 6176 bool calc_cond) const
6177 { 6177 {
6178 ComplexMatrix retval; 6178 ComplexMatrix retval;
6188 { 6188 {
6189 // Print spparms("spumoni") info if requested 6189 // Print spparms("spumoni") info if requested
6190 volatile int typ = mattype.type (); 6190 volatile int typ = mattype.type ();
6191 mattype.info (); 6191 mattype.info ();
6192 6192
6193 if (typ == SparseType::Hermitian) 6193 if (typ == MatrixType::Hermitian)
6194 { 6194 {
6195 #ifdef HAVE_CHOLMOD 6195 #ifdef HAVE_CHOLMOD
6196 cholmod_common Common; 6196 cholmod_common Common;
6197 cholmod_common *cm = &Common; 6197 cholmod_common *cm = &Common;
6198 6198
6290 6290
6291 if (rcond == 0.0) 6291 if (rcond == 0.0)
6292 { 6292 {
6293 // Either its indefinite or singular. Try UMFPACK 6293 // Either its indefinite or singular. Try UMFPACK
6294 mattype.mark_as_unsymmetric (); 6294 mattype.mark_as_unsymmetric ();
6295 typ = SparseType::Full; 6295 typ = MatrixType::Full;
6296 } 6296 }
6297 else 6297 else
6298 { 6298 {
6299 volatile double rcond_plus_one = rcond + 1.0; 6299 volatile double rcond_plus_one = rcond + 1.0;
6300 6300
6338 #else 6338 #else
6339 (*current_liboctave_warning_handler) 6339 (*current_liboctave_warning_handler)
6340 ("CHOLMOD not installed"); 6340 ("CHOLMOD not installed");
6341 6341
6342 mattype.mark_as_unsymmetric (); 6342 mattype.mark_as_unsymmetric ();
6343 typ = SparseType::Full; 6343 typ = MatrixType::Full;
6344 #endif 6344 #endif
6345 } 6345 }
6346 6346
6347 if (typ == SparseType::Full) 6347 if (typ == MatrixType::Full)
6348 { 6348 {
6349 #ifdef HAVE_UMFPACK 6349 #ifdef HAVE_UMFPACK
6350 Matrix Control, Info; 6350 Matrix Control, Info;
6351 void *Numeric = factorize (err, rcond, Control, Info, 6351 void *Numeric = factorize (err, rcond, Control, Info,
6352 sing_handler, calc_cond); 6352 sing_handler, calc_cond);
6399 6399
6400 #else 6400 #else
6401 (*current_liboctave_error_handler) ("UMFPACK not installed"); 6401 (*current_liboctave_error_handler) ("UMFPACK not installed");
6402 #endif 6402 #endif
6403 } 6403 }
6404 else if (typ != SparseType::Hermitian) 6404 else if (typ != MatrixType::Hermitian)
6405 (*current_liboctave_error_handler) ("incorrect matrix type"); 6405 (*current_liboctave_error_handler) ("incorrect matrix type");
6406 } 6406 }
6407 6407
6408 return retval; 6408 return retval;
6409 } 6409 }
6410 6410
6411 SparseComplexMatrix 6411 SparseComplexMatrix
6412 SparseComplexMatrix::fsolve (SparseType &mattype, const SparseComplexMatrix& b, 6412 SparseComplexMatrix::fsolve (MatrixType &mattype, const SparseComplexMatrix& b,
6413 octave_idx_type& err, double& rcond, 6413 octave_idx_type& err, double& rcond,
6414 solve_singularity_handler sing_handler, 6414 solve_singularity_handler sing_handler,
6415 bool calc_cond) const 6415 bool calc_cond) const
6416 { 6416 {
6417 SparseComplexMatrix retval; 6417 SparseComplexMatrix retval;
6427 { 6427 {
6428 // Print spparms("spumoni") info if requested 6428 // Print spparms("spumoni") info if requested
6429 volatile int typ = mattype.type (); 6429 volatile int typ = mattype.type ();
6430 mattype.info (); 6430 mattype.info ();
6431 6431
6432 if (typ == SparseType::Hermitian) 6432 if (typ == MatrixType::Hermitian)
6433 { 6433 {
6434 #ifdef HAVE_CHOLMOD 6434 #ifdef HAVE_CHOLMOD
6435 cholmod_common Common; 6435 cholmod_common Common;
6436 cholmod_common *cm = &Common; 6436 cholmod_common *cm = &Common;
6437 6437
6539 6539
6540 if (rcond == 0.0) 6540 if (rcond == 0.0)
6541 { 6541 {
6542 // Either its indefinite or singular. Try UMFPACK 6542 // Either its indefinite or singular. Try UMFPACK
6543 mattype.mark_as_unsymmetric (); 6543 mattype.mark_as_unsymmetric ();
6544 typ = SparseType::Full; 6544 typ = MatrixType::Full;
6545 } 6545 }
6546 else 6546 else
6547 { 6547 {
6548 volatile double rcond_plus_one = rcond + 1.0; 6548 volatile double rcond_plus_one = rcond + 1.0;
6549 6549
6593 #else 6593 #else
6594 (*current_liboctave_warning_handler) 6594 (*current_liboctave_warning_handler)
6595 ("CHOLMOD not installed"); 6595 ("CHOLMOD not installed");
6596 6596
6597 mattype.mark_as_unsymmetric (); 6597 mattype.mark_as_unsymmetric ();
6598 typ = SparseType::Full; 6598 typ = MatrixType::Full;
6599 #endif 6599 #endif
6600 } 6600 }
6601 6601
6602 if (typ == SparseType::Full) 6602 if (typ == MatrixType::Full)
6603 { 6603 {
6604 #ifdef HAVE_UMFPACK 6604 #ifdef HAVE_UMFPACK
6605 Matrix Control, Info; 6605 Matrix Control, Info;
6606 void *Numeric = factorize (err, rcond, Control, Info, 6606 void *Numeric = factorize (err, rcond, Control, Info,
6607 sing_handler, calc_cond); 6607 sing_handler, calc_cond);
6702 6702
6703 #else 6703 #else
6704 (*current_liboctave_error_handler) ("UMFPACK not installed"); 6704 (*current_liboctave_error_handler) ("UMFPACK not installed");
6705 #endif 6705 #endif
6706 } 6706 }
6707 else if (typ != SparseType::Hermitian) 6707 else if (typ != MatrixType::Hermitian)
6708 (*current_liboctave_error_handler) ("incorrect matrix type"); 6708 (*current_liboctave_error_handler) ("incorrect matrix type");
6709 } 6709 }
6710 6710
6711 return retval; 6711 return retval;
6712 } 6712 }
6713 6713
6714 ComplexMatrix 6714 ComplexMatrix
6715 SparseComplexMatrix::solve (SparseType &mattype, const Matrix& b) const 6715 SparseComplexMatrix::solve (MatrixType &mattype, const Matrix& b) const
6716 { 6716 {
6717 octave_idx_type info; 6717 octave_idx_type info;
6718 double rcond; 6718 double rcond;
6719 return solve (mattype, b, info, rcond, 0); 6719 return solve (mattype, b, info, rcond, 0);
6720 } 6720 }
6721 6721
6722 ComplexMatrix 6722 ComplexMatrix
6723 SparseComplexMatrix::solve (SparseType &mattype, const Matrix& b, 6723 SparseComplexMatrix::solve (MatrixType &mattype, const Matrix& b,
6724 octave_idx_type& info) const 6724 octave_idx_type& info) const
6725 { 6725 {
6726 double rcond; 6726 double rcond;
6727 return solve (mattype, b, info, rcond, 0); 6727 return solve (mattype, b, info, rcond, 0);
6728 } 6728 }
6729 6729
6730 ComplexMatrix 6730 ComplexMatrix
6731 SparseComplexMatrix::solve (SparseType &mattype, const Matrix& b, 6731 SparseComplexMatrix::solve (MatrixType &mattype, const Matrix& b,
6732 octave_idx_type& info, double& rcond) const 6732 octave_idx_type& info, double& rcond) const
6733 { 6733 {
6734 return solve (mattype, b, info, rcond, 0); 6734 return solve (mattype, b, info, rcond, 0);
6735 } 6735 }
6736 6736
6737 ComplexMatrix 6737 ComplexMatrix
6738 SparseComplexMatrix::solve (SparseType &mattype, const Matrix& b, 6738 SparseComplexMatrix::solve (MatrixType &mattype, const Matrix& b,
6739 octave_idx_type& err, double& rcond, 6739 octave_idx_type& err, double& rcond,
6740 solve_singularity_handler sing_handler, 6740 solve_singularity_handler sing_handler,
6741 bool singular_fallback) const 6741 bool singular_fallback) const
6742 { 6742 {
6743 ComplexMatrix retval; 6743 ComplexMatrix retval;
6744 int typ = mattype.type (false); 6744 int typ = mattype.type (false);
6745 6745
6746 if (typ == SparseType::Unknown) 6746 if (typ == MatrixType::Unknown)
6747 typ = mattype.type (*this); 6747 typ = mattype.type (*this);
6748 6748
6749 if (typ == SparseType::Diagonal || typ == SparseType::Permuted_Diagonal) 6749 if (typ == MatrixType::Diagonal || typ == MatrixType::Permuted_Diagonal)
6750 retval = dsolve (mattype, b, err, rcond, sing_handler, false); 6750 retval = dsolve (mattype, b, err, rcond, sing_handler, false);
6751 else if (typ == SparseType::Upper || typ == SparseType::Permuted_Upper) 6751 else if (typ == MatrixType::Upper || typ == MatrixType::Permuted_Upper)
6752 retval = utsolve (mattype, b, err, rcond, sing_handler, false); 6752 retval = utsolve (mattype, b, err, rcond, sing_handler, false);
6753 else if (typ == SparseType::Lower || typ == SparseType::Permuted_Lower) 6753 else if (typ == MatrixType::Lower || typ == MatrixType::Permuted_Lower)
6754 retval = ltsolve (mattype, b, err, rcond, sing_handler, false); 6754 retval = ltsolve (mattype, b, err, rcond, sing_handler, false);
6755 else if (typ == SparseType::Banded || typ == SparseType::Banded_Hermitian) 6755 else if (typ == MatrixType::Banded || typ == MatrixType::Banded_Hermitian)
6756 retval = bsolve (mattype, b, err, rcond, sing_handler, false); 6756 retval = bsolve (mattype, b, err, rcond, sing_handler, false);
6757 else if (typ == SparseType::Tridiagonal || 6757 else if (typ == MatrixType::Tridiagonal ||
6758 typ == SparseType::Tridiagonal_Hermitian) 6758 typ == MatrixType::Tridiagonal_Hermitian)
6759 retval = trisolve (mattype, b, err, rcond, sing_handler, false); 6759 retval = trisolve (mattype, b, err, rcond, sing_handler, false);
6760 else if (typ == SparseType::Full || typ == SparseType::Hermitian) 6760 else if (typ == MatrixType::Full || typ == MatrixType::Hermitian)
6761 retval = fsolve (mattype, b, err, rcond, sing_handler, true); 6761 retval = fsolve (mattype, b, err, rcond, sing_handler, true);
6762 else if (typ != SparseType::Rectangular) 6762 else if (typ != MatrixType::Rectangular)
6763 { 6763 {
6764 (*current_liboctave_error_handler) ("unknown matrix type"); 6764 (*current_liboctave_error_handler) ("unknown matrix type");
6765 return ComplexMatrix (); 6765 return ComplexMatrix ();
6766 } 6766 }
6767 6767
6768 if (singular_fallback && mattype.type(false) == SparseType::Rectangular) 6768 if (singular_fallback && mattype.type(false) == MatrixType::Rectangular)
6769 { 6769 {
6770 rcond = 1.; 6770 rcond = 1.;
6771 #ifdef USE_QRSOLVE 6771 #ifdef USE_QRSOLVE
6772 retval = qrsolve (*this, b, err); 6772 retval = qrsolve (*this, b, err);
6773 #else 6773 #else
6778 6778
6779 return retval; 6779 return retval;
6780 } 6780 }
6781 6781
6782 SparseComplexMatrix 6782 SparseComplexMatrix
6783 SparseComplexMatrix::solve (SparseType &mattype, const SparseMatrix& b) const 6783 SparseComplexMatrix::solve (MatrixType &mattype, const SparseMatrix& b) const
6784 { 6784 {
6785 octave_idx_type info; 6785 octave_idx_type info;
6786 double rcond; 6786 double rcond;
6787 return solve (mattype, b, info, rcond, 0); 6787 return solve (mattype, b, info, rcond, 0);
6788 } 6788 }
6789 6789
6790 SparseComplexMatrix 6790 SparseComplexMatrix
6791 SparseComplexMatrix::solve (SparseType &mattype, const SparseMatrix& b, 6791 SparseComplexMatrix::solve (MatrixType &mattype, const SparseMatrix& b,
6792 octave_idx_type& info) const 6792 octave_idx_type& info) const
6793 { 6793 {
6794 double rcond; 6794 double rcond;
6795 return solve (mattype, b, info, rcond, 0); 6795 return solve (mattype, b, info, rcond, 0);
6796 } 6796 }
6797 6797
6798 SparseComplexMatrix 6798 SparseComplexMatrix
6799 SparseComplexMatrix::solve (SparseType &mattype, const SparseMatrix& b, 6799 SparseComplexMatrix::solve (MatrixType &mattype, const SparseMatrix& b,
6800 octave_idx_type& info, double& rcond) const 6800 octave_idx_type& info, double& rcond) const
6801 { 6801 {
6802 return solve (mattype, b, info, rcond, 0); 6802 return solve (mattype, b, info, rcond, 0);
6803 } 6803 }
6804 6804
6805 SparseComplexMatrix 6805 SparseComplexMatrix
6806 SparseComplexMatrix::solve (SparseType &mattype, const SparseMatrix& b, 6806 SparseComplexMatrix::solve (MatrixType &mattype, const SparseMatrix& b,
6807 octave_idx_type& err, double& rcond, 6807 octave_idx_type& err, double& rcond,
6808 solve_singularity_handler sing_handler, 6808 solve_singularity_handler sing_handler,
6809 bool singular_fallback) const 6809 bool singular_fallback) const
6810 { 6810 {
6811 SparseComplexMatrix retval; 6811 SparseComplexMatrix retval;
6812 int typ = mattype.type (false); 6812 int typ = mattype.type (false);
6813 6813
6814 if (typ == SparseType::Unknown) 6814 if (typ == MatrixType::Unknown)
6815 typ = mattype.type (*this); 6815 typ = mattype.type (*this);
6816 6816
6817 if (typ == SparseType::Diagonal || typ == SparseType::Permuted_Diagonal) 6817 if (typ == MatrixType::Diagonal || typ == MatrixType::Permuted_Diagonal)
6818 retval = dsolve (mattype, b, err, rcond, sing_handler, false); 6818 retval = dsolve (mattype, b, err, rcond, sing_handler, false);
6819 else if (typ == SparseType::Upper || typ == SparseType::Permuted_Upper) 6819 else if (typ == MatrixType::Upper || typ == MatrixType::Permuted_Upper)
6820 retval = utsolve (mattype, b, err, rcond, sing_handler, false); 6820 retval = utsolve (mattype, b, err, rcond, sing_handler, false);
6821 else if (typ == SparseType::Lower || typ == SparseType::Permuted_Lower) 6821 else if (typ == MatrixType::Lower || typ == MatrixType::Permuted_Lower)
6822 retval = ltsolve (mattype, b, err, rcond, sing_handler, false); 6822 retval = ltsolve (mattype, b, err, rcond, sing_handler, false);
6823 else if (typ == SparseType::Banded || typ == SparseType::Banded_Hermitian) 6823 else if (typ == MatrixType::Banded || typ == MatrixType::Banded_Hermitian)
6824 retval = bsolve (mattype, b, err, rcond, sing_handler, false); 6824 retval = bsolve (mattype, b, err, rcond, sing_handler, false);
6825 else if (typ == SparseType::Tridiagonal || 6825 else if (typ == MatrixType::Tridiagonal ||
6826 typ == SparseType::Tridiagonal_Hermitian) 6826 typ == MatrixType::Tridiagonal_Hermitian)
6827 retval = trisolve (mattype, b, err, rcond, sing_handler, false); 6827 retval = trisolve (mattype, b, err, rcond, sing_handler, false);
6828 else if (typ == SparseType::Full || typ == SparseType::Hermitian) 6828 else if (typ == MatrixType::Full || typ == MatrixType::Hermitian)
6829 retval = fsolve (mattype, b, err, rcond, sing_handler, true); 6829 retval = fsolve (mattype, b, err, rcond, sing_handler, true);
6830 else if (typ != SparseType::Rectangular) 6830 else if (typ != MatrixType::Rectangular)
6831 { 6831 {
6832 (*current_liboctave_error_handler) ("unknown matrix type"); 6832 (*current_liboctave_error_handler) ("unknown matrix type");
6833 return SparseComplexMatrix (); 6833 return SparseComplexMatrix ();
6834 } 6834 }
6835 6835
6836 if (singular_fallback && mattype.type(false) == SparseType::Rectangular) 6836 if (singular_fallback && mattype.type(false) == MatrixType::Rectangular)
6837 { 6837 {
6838 rcond = 1.; 6838 rcond = 1.;
6839 #ifdef USE_QRSOLVE 6839 #ifdef USE_QRSOLVE
6840 retval = qrsolve (*this, b, err); 6840 retval = qrsolve (*this, b, err);
6841 #else 6841 #else
6846 6846
6847 return retval; 6847 return retval;
6848 } 6848 }
6849 6849
6850 ComplexMatrix 6850 ComplexMatrix
6851 SparseComplexMatrix::solve (SparseType &mattype, const ComplexMatrix& b) const 6851 SparseComplexMatrix::solve (MatrixType &mattype, const ComplexMatrix& b) const
6852 { 6852 {
6853 octave_idx_type info; 6853 octave_idx_type info;
6854 double rcond; 6854 double rcond;
6855 return solve (mattype, b, info, rcond, 0); 6855 return solve (mattype, b, info, rcond, 0);
6856 } 6856 }
6857 6857
6858 ComplexMatrix 6858 ComplexMatrix
6859 SparseComplexMatrix::solve (SparseType &mattype, const ComplexMatrix& b, 6859 SparseComplexMatrix::solve (MatrixType &mattype, const ComplexMatrix& b,
6860 octave_idx_type& info) const 6860 octave_idx_type& info) const
6861 { 6861 {
6862 double rcond; 6862 double rcond;
6863 return solve (mattype, b, info, rcond, 0); 6863 return solve (mattype, b, info, rcond, 0);
6864 } 6864 }
6865 6865
6866 ComplexMatrix 6866 ComplexMatrix
6867 SparseComplexMatrix::solve (SparseType &mattype, const ComplexMatrix& b, 6867 SparseComplexMatrix::solve (MatrixType &mattype, const ComplexMatrix& b,
6868 octave_idx_type& info, double& rcond) const 6868 octave_idx_type& info, double& rcond) const
6869 { 6869 {
6870 return solve (mattype, b, info, rcond, 0); 6870 return solve (mattype, b, info, rcond, 0);
6871 } 6871 }
6872 6872
6873 ComplexMatrix 6873 ComplexMatrix
6874 SparseComplexMatrix::solve (SparseType &mattype, const ComplexMatrix& b, 6874 SparseComplexMatrix::solve (MatrixType &mattype, const ComplexMatrix& b,
6875 octave_idx_type& err, double& rcond, 6875 octave_idx_type& err, double& rcond,
6876 solve_singularity_handler sing_handler, 6876 solve_singularity_handler sing_handler,
6877 bool singular_fallback) const 6877 bool singular_fallback) const
6878 { 6878 {
6879 ComplexMatrix retval; 6879 ComplexMatrix retval;
6880 int typ = mattype.type (false); 6880 int typ = mattype.type (false);
6881 6881
6882 if (typ == SparseType::Unknown) 6882 if (typ == MatrixType::Unknown)
6883 typ = mattype.type (*this); 6883 typ = mattype.type (*this);
6884 6884
6885 if (typ == SparseType::Diagonal || typ == SparseType::Permuted_Diagonal) 6885 if (typ == MatrixType::Diagonal || typ == MatrixType::Permuted_Diagonal)
6886 retval = dsolve (mattype, b, err, rcond, sing_handler, false); 6886 retval = dsolve (mattype, b, err, rcond, sing_handler, false);
6887 else if (typ == SparseType::Upper || typ == SparseType::Permuted_Upper) 6887 else if (typ == MatrixType::Upper || typ == MatrixType::Permuted_Upper)
6888 retval = utsolve (mattype, b, err, rcond, sing_handler, false); 6888 retval = utsolve (mattype, b, err, rcond, sing_handler, false);
6889 else if (typ == SparseType::Lower || typ == SparseType::Permuted_Lower) 6889 else if (typ == MatrixType::Lower || typ == MatrixType::Permuted_Lower)
6890 retval = ltsolve (mattype, b, err, rcond, sing_handler, false); 6890 retval = ltsolve (mattype, b, err, rcond, sing_handler, false);
6891 else if (typ == SparseType::Banded || typ == SparseType::Banded_Hermitian) 6891 else if (typ == MatrixType::Banded || typ == MatrixType::Banded_Hermitian)
6892 retval = bsolve (mattype, b, err, rcond, sing_handler, false); 6892 retval = bsolve (mattype, b, err, rcond, sing_handler, false);
6893 else if (typ == SparseType::Tridiagonal || 6893 else if (typ == MatrixType::Tridiagonal ||
6894 typ == SparseType::Tridiagonal_Hermitian) 6894 typ == MatrixType::Tridiagonal_Hermitian)
6895 retval = trisolve (mattype, b, err, rcond, sing_handler, false); 6895 retval = trisolve (mattype, b, err, rcond, sing_handler, false);
6896 else if (typ == SparseType::Full || typ == SparseType::Hermitian) 6896 else if (typ == MatrixType::Full || typ == MatrixType::Hermitian)
6897 retval = fsolve (mattype, b, err, rcond, sing_handler, true); 6897 retval = fsolve (mattype, b, err, rcond, sing_handler, true);
6898 else if (typ != SparseType::Rectangular) 6898 else if (typ != MatrixType::Rectangular)
6899 { 6899 {
6900 (*current_liboctave_error_handler) ("unknown matrix type"); 6900 (*current_liboctave_error_handler) ("unknown matrix type");
6901 return ComplexMatrix (); 6901 return ComplexMatrix ();
6902 } 6902 }
6903 6903
6904 if (singular_fallback && mattype.type(false) == SparseType::Rectangular) 6904 if (singular_fallback && mattype.type(false) == MatrixType::Rectangular)
6905 { 6905 {
6906 rcond = 1.; 6906 rcond = 1.;
6907 #ifdef USE_QRSOLVE 6907 #ifdef USE_QRSOLVE
6908 retval = qrsolve (*this, b, err); 6908 retval = qrsolve (*this, b, err);
6909 #else 6909 #else
6914 6914
6915 return retval; 6915 return retval;
6916 } 6916 }
6917 6917
6918 SparseComplexMatrix 6918 SparseComplexMatrix
6919 SparseComplexMatrix::solve (SparseType &mattype, 6919 SparseComplexMatrix::solve (MatrixType &mattype,
6920 const SparseComplexMatrix& b) const 6920 const SparseComplexMatrix& b) const
6921 { 6921 {
6922 octave_idx_type info; 6922 octave_idx_type info;
6923 double rcond; 6923 double rcond;
6924 return solve (mattype, b, info, rcond, 0); 6924 return solve (mattype, b, info, rcond, 0);
6925 } 6925 }
6926 6926
6927 SparseComplexMatrix 6927 SparseComplexMatrix
6928 SparseComplexMatrix::solve (SparseType &mattype, const SparseComplexMatrix& b, 6928 SparseComplexMatrix::solve (MatrixType &mattype, const SparseComplexMatrix& b,
6929 octave_idx_type& info) const 6929 octave_idx_type& info) const
6930 { 6930 {
6931 double rcond; 6931 double rcond;
6932 return solve (mattype, b, info, rcond, 0); 6932 return solve (mattype, b, info, rcond, 0);
6933 } 6933 }
6934 6934
6935 SparseComplexMatrix 6935 SparseComplexMatrix
6936 SparseComplexMatrix::solve (SparseType &mattype, const SparseComplexMatrix& b, 6936 SparseComplexMatrix::solve (MatrixType &mattype, const SparseComplexMatrix& b,
6937 octave_idx_type& info, double& rcond) const 6937 octave_idx_type& info, double& rcond) const
6938 { 6938 {
6939 return solve (mattype, b, info, rcond, 0); 6939 return solve (mattype, b, info, rcond, 0);
6940 } 6940 }
6941 6941
6942 SparseComplexMatrix 6942 SparseComplexMatrix
6943 SparseComplexMatrix::solve (SparseType &mattype, const SparseComplexMatrix& b, 6943 SparseComplexMatrix::solve (MatrixType &mattype, const SparseComplexMatrix& b,
6944 octave_idx_type& err, double& rcond, 6944 octave_idx_type& err, double& rcond,
6945 solve_singularity_handler sing_handler, 6945 solve_singularity_handler sing_handler,
6946 bool singular_fallback) const 6946 bool singular_fallback) const
6947 { 6947 {
6948 SparseComplexMatrix retval; 6948 SparseComplexMatrix retval;
6949 int typ = mattype.type (false); 6949 int typ = mattype.type (false);
6950 6950
6951 if (typ == SparseType::Unknown) 6951 if (typ == MatrixType::Unknown)
6952 typ = mattype.type (*this); 6952 typ = mattype.type (*this);
6953 6953
6954 if (typ == SparseType::Diagonal || typ == SparseType::Permuted_Diagonal) 6954 if (typ == MatrixType::Diagonal || typ == MatrixType::Permuted_Diagonal)
6955 retval = dsolve (mattype, b, err, rcond, sing_handler, false); 6955 retval = dsolve (mattype, b, err, rcond, sing_handler, false);
6956 else if (typ == SparseType::Upper || typ == SparseType::Permuted_Upper) 6956 else if (typ == MatrixType::Upper || typ == MatrixType::Permuted_Upper)
6957 retval = utsolve (mattype, b, err, rcond, sing_handler, false); 6957 retval = utsolve (mattype, b, err, rcond, sing_handler, false);
6958 else if (typ == SparseType::Lower || typ == SparseType::Permuted_Lower) 6958 else if (typ == MatrixType::Lower || typ == MatrixType::Permuted_Lower)
6959 retval = ltsolve (mattype, b, err, rcond, sing_handler, false); 6959 retval = ltsolve (mattype, b, err, rcond, sing_handler, false);
6960 else if (typ == SparseType::Banded || typ == SparseType::Banded_Hermitian) 6960 else if (typ == MatrixType::Banded || typ == MatrixType::Banded_Hermitian)
6961 retval = bsolve (mattype, b, err, rcond, sing_handler, false); 6961 retval = bsolve (mattype, b, err, rcond, sing_handler, false);
6962 else if (typ == SparseType::Tridiagonal || 6962 else if (typ == MatrixType::Tridiagonal ||
6963 typ == SparseType::Tridiagonal_Hermitian) 6963 typ == MatrixType::Tridiagonal_Hermitian)
6964 retval = trisolve (mattype, b, err, rcond, sing_handler, false); 6964 retval = trisolve (mattype, b, err, rcond, sing_handler, false);
6965 else if (typ == SparseType::Full || typ == SparseType::Hermitian) 6965 else if (typ == MatrixType::Full || typ == MatrixType::Hermitian)
6966 retval = fsolve (mattype, b, err, rcond, sing_handler, true); 6966 retval = fsolve (mattype, b, err, rcond, sing_handler, true);
6967 else if (typ != SparseType::Rectangular) 6967 else if (typ != MatrixType::Rectangular)
6968 { 6968 {
6969 (*current_liboctave_error_handler) ("unknown matrix type"); 6969 (*current_liboctave_error_handler) ("unknown matrix type");
6970 return SparseComplexMatrix (); 6970 return SparseComplexMatrix ();
6971 } 6971 }
6972 6972
6973 if (singular_fallback && mattype.type(false) == SparseType::Rectangular) 6973 if (singular_fallback && mattype.type(false) == MatrixType::Rectangular)
6974 { 6974 {
6975 rcond = 1.; 6975 rcond = 1.;
6976 #ifdef USE_QRSOLVE 6976 #ifdef USE_QRSOLVE
6977 retval = qrsolve (*this, b, err); 6977 retval = qrsolve (*this, b, err);
6978 #else 6978 #else
6983 6983
6984 return retval; 6984 return retval;
6985 } 6985 }
6986 6986
6987 ComplexColumnVector 6987 ComplexColumnVector
6988 SparseComplexMatrix::solve (SparseType &mattype, const ColumnVector& b) const 6988 SparseComplexMatrix::solve (MatrixType &mattype, const ColumnVector& b) const
6989 { 6989 {
6990 octave_idx_type info; double rcond; 6990 octave_idx_type info; double rcond;
6991 return solve (mattype, b, info, rcond); 6991 return solve (mattype, b, info, rcond);
6992 } 6992 }
6993 6993
6994 ComplexColumnVector 6994 ComplexColumnVector
6995 SparseComplexMatrix::solve (SparseType &mattype, const ColumnVector& b, 6995 SparseComplexMatrix::solve (MatrixType &mattype, const ColumnVector& b,
6996 octave_idx_type& info) const 6996 octave_idx_type& info) const
6997 { 6997 {
6998 double rcond; 6998 double rcond;
6999 return solve (mattype, b, info, rcond); 6999 return solve (mattype, b, info, rcond);
7000 } 7000 }
7001 7001
7002 ComplexColumnVector 7002 ComplexColumnVector
7003 SparseComplexMatrix::solve (SparseType &mattype, const ColumnVector& b, 7003 SparseComplexMatrix::solve (MatrixType &mattype, const ColumnVector& b,
7004 octave_idx_type& info, double& rcond) const 7004 octave_idx_type& info, double& rcond) const
7005 { 7005 {
7006 return solve (mattype, b, info, rcond, 0); 7006 return solve (mattype, b, info, rcond, 0);
7007 } 7007 }
7008 7008
7009 ComplexColumnVector 7009 ComplexColumnVector
7010 SparseComplexMatrix::solve (SparseType &mattype, const ColumnVector& b, 7010 SparseComplexMatrix::solve (MatrixType &mattype, const ColumnVector& b,
7011 octave_idx_type& info, double& rcond, 7011 octave_idx_type& info, double& rcond,
7012 solve_singularity_handler sing_handler) const 7012 solve_singularity_handler sing_handler) const
7013 { 7013 {
7014 Matrix tmp (b); 7014 Matrix tmp (b);
7015 return solve (mattype, tmp, info, rcond, sing_handler).column (static_cast<octave_idx_type> (0)); 7015 return solve (mattype, tmp, info, rcond, sing_handler).column (static_cast<octave_idx_type> (0));
7016 } 7016 }
7017 7017
7018 ComplexColumnVector 7018 ComplexColumnVector
7019 SparseComplexMatrix::solve (SparseType &mattype, 7019 SparseComplexMatrix::solve (MatrixType &mattype,
7020 const ComplexColumnVector& b) const 7020 const ComplexColumnVector& b) const
7021 { 7021 {
7022 octave_idx_type info; 7022 octave_idx_type info;
7023 double rcond; 7023 double rcond;
7024 return solve (mattype, b, info, rcond, 0); 7024 return solve (mattype, b, info, rcond, 0);
7025 } 7025 }
7026 7026
7027 ComplexColumnVector 7027 ComplexColumnVector
7028 SparseComplexMatrix::solve (SparseType &mattype, const ComplexColumnVector& b, 7028 SparseComplexMatrix::solve (MatrixType &mattype, const ComplexColumnVector& b,
7029 octave_idx_type& info) const 7029 octave_idx_type& info) const
7030 { 7030 {
7031 double rcond; 7031 double rcond;
7032 return solve (mattype, b, info, rcond, 0); 7032 return solve (mattype, b, info, rcond, 0);
7033 } 7033 }
7034 7034
7035 ComplexColumnVector 7035 ComplexColumnVector
7036 SparseComplexMatrix::solve (SparseType &mattype, const ComplexColumnVector& b, 7036 SparseComplexMatrix::solve (MatrixType &mattype, const ComplexColumnVector& b,
7037 octave_idx_type& info, double& rcond) const 7037 octave_idx_type& info, double& rcond) const
7038 { 7038 {
7039 return solve (mattype, b, info, rcond, 0); 7039 return solve (mattype, b, info, rcond, 0);
7040 } 7040 }
7041 7041
7042 ComplexColumnVector 7042 ComplexColumnVector
7043 SparseComplexMatrix::solve (SparseType &mattype, const ComplexColumnVector& b, 7043 SparseComplexMatrix::solve (MatrixType &mattype, const ComplexColumnVector& b,
7044 octave_idx_type& info, double& rcond, 7044 octave_idx_type& info, double& rcond,
7045 solve_singularity_handler sing_handler) const 7045 solve_singularity_handler sing_handler) const
7046 { 7046 {
7047 ComplexMatrix tmp (b); 7047 ComplexMatrix tmp (b);
7048 return solve (mattype, tmp, info, rcond, sing_handler).column (static_cast<octave_idx_type> (0)); 7048 return solve (mattype, tmp, info, rcond, sing_handler).column (static_cast<octave_idx_type> (0));
7073 ComplexMatrix 7073 ComplexMatrix
7074 SparseComplexMatrix::solve (const Matrix& b, octave_idx_type& err, 7074 SparseComplexMatrix::solve (const Matrix& b, octave_idx_type& err,
7075 double& rcond, 7075 double& rcond,
7076 solve_singularity_handler sing_handler) const 7076 solve_singularity_handler sing_handler) const
7077 { 7077 {
7078 SparseType mattype (*this); 7078 MatrixType mattype (*this);
7079 return solve (mattype, b, err, rcond, sing_handler); 7079 return solve (mattype, b, err, rcond, sing_handler);
7080 } 7080 }
7081 7081
7082 SparseComplexMatrix 7082 SparseComplexMatrix
7083 SparseComplexMatrix::solve (const SparseMatrix& b) const 7083 SparseComplexMatrix::solve (const SparseMatrix& b) const
7105 SparseComplexMatrix 7105 SparseComplexMatrix
7106 SparseComplexMatrix::solve (const SparseMatrix& b, 7106 SparseComplexMatrix::solve (const SparseMatrix& b,
7107 octave_idx_type& err, double& rcond, 7107 octave_idx_type& err, double& rcond,
7108 solve_singularity_handler sing_handler) const 7108 solve_singularity_handler sing_handler) const
7109 { 7109 {
7110 SparseType mattype (*this); 7110 MatrixType mattype (*this);
7111 return solve (mattype, b, err, rcond, sing_handler); 7111 return solve (mattype, b, err, rcond, sing_handler);
7112 } 7112 }
7113 7113
7114 ComplexMatrix 7114 ComplexMatrix
7115 SparseComplexMatrix::solve (const ComplexMatrix& b, 7115 SparseComplexMatrix::solve (const ComplexMatrix& b,
7129 ComplexMatrix 7129 ComplexMatrix
7130 SparseComplexMatrix::solve (const ComplexMatrix& b, 7130 SparseComplexMatrix::solve (const ComplexMatrix& b,
7131 octave_idx_type& err, double& rcond, 7131 octave_idx_type& err, double& rcond,
7132 solve_singularity_handler sing_handler) const 7132 solve_singularity_handler sing_handler) const
7133 { 7133 {
7134 SparseType mattype (*this); 7134 MatrixType mattype (*this);
7135 return solve (mattype, b, err, rcond, sing_handler); 7135 return solve (mattype, b, err, rcond, sing_handler);
7136 } 7136 }
7137 7137
7138 SparseComplexMatrix 7138 SparseComplexMatrix
7139 SparseComplexMatrix::solve (const SparseComplexMatrix& b) const 7139 SparseComplexMatrix::solve (const SparseComplexMatrix& b) const
7161 SparseComplexMatrix 7161 SparseComplexMatrix
7162 SparseComplexMatrix::solve (const SparseComplexMatrix& b, 7162 SparseComplexMatrix::solve (const SparseComplexMatrix& b,
7163 octave_idx_type& err, double& rcond, 7163 octave_idx_type& err, double& rcond,
7164 solve_singularity_handler sing_handler) const 7164 solve_singularity_handler sing_handler) const
7165 { 7165 {
7166 SparseType mattype (*this); 7166 MatrixType mattype (*this);
7167 return solve (mattype, b, err, rcond, sing_handler); 7167 return solve (mattype, b, err, rcond, sing_handler);
7168 } 7168 }
7169 7169
7170 ComplexColumnVector 7170 ComplexColumnVector
7171 SparseComplexMatrix::solve (const ColumnVector& b) const 7171 SparseComplexMatrix::solve (const ColumnVector& b) const