26 #if defined (HAVE_CONFIG_H)
69 #if ! defined (USE_QRSOLVE)
121 if (nr != nr_a || nc != nc_a || nz != nz_a)
138 return !(*
this == a);
147 if (nr == nc && nr > 0)
189 return max (dummy_idx, dim);
200 if (dim >= dv.
ndims ())
213 if (nr == 0 || nc == 0 || dim >= dv.
ndims ())
224 if (
ridx (i) != idx_j)
243 double abs_tmp = std::abs (tmp);
261 result.
xcidx (0) = 0;
267 result.
xdata (ii) = tmp;
268 result.
xridx (ii++) = 0;
270 result.
xcidx (j+1) = ii;
277 if (nr == 0 || nc == 0 || dim >= dv.
ndims ())
287 if (found[
ridx (i)] == -j)
288 found[
ridx (i)] = -j - 1;
291 if (found[i] > -nc && found[i] < 0)
292 idx_arg.
elem (i) = -found[i];
304 else if (ix == -1 || std::abs (tmp) > std::abs (
elem (ir, ix)))
305 idx_arg.
elem (ir) = j;
311 if (idx_arg.
elem (j) == -1 ||
elem (j, idx_arg.
elem (j)) != 0.)
317 result.
xcidx (0) = 0;
318 result.
xcidx (1) = nel;
321 if (idx_arg(j) == -1)
324 result.
xdata (ii) = Complex_NaN_result;
325 result.
xridx (ii++) = j;
332 result.
xdata (ii) = tmp;
333 result.
xridx (ii++) = j;
346 return min (dummy_idx, dim);
357 if (dim >= dv.
ndims ())
370 if (nr == 0 || nc == 0 || dim >= dv.
ndims ())
381 if (
ridx (i) != idx_j)
400 double abs_tmp = std::abs (tmp);
418 result.
xcidx (0) = 0;
424 result.
xdata (ii) = tmp;
425 result.
xridx (ii++) = 0;
427 result.
xcidx (j+1) = ii;
434 if (nr == 0 || nc == 0 || dim >= dv.
ndims ())
444 if (found[
ridx (i)] == -j)
445 found[
ridx (i)] = -j - 1;
448 if (found[i] > -nc && found[i] < 0)
449 idx_arg.
elem (i) = -found[i];
461 else if (ix == -1 || std::abs (tmp) < std::abs (
elem (ir, ix)))
462 idx_arg.
elem (ir) = j;
468 if (idx_arg.
elem (j) == -1 ||
elem (j, idx_arg.
elem (j)) != 0.)
474 result.
xcidx (0) = 0;
475 result.
xcidx (1) = nel;
478 if (idx_arg(j) == -1)
481 result.
xdata (ii) = Complex_NaN_result;
482 result.
xridx (ii++) = j;
489 result.
xdata (ii) = tmp;
490 result.
xridx (ii++) = j;
525 retval(j) =
data (k);
568 return insert (tmp , indx);
584 if (rb.
rows () > 0 && rb.
cols () > 0)
594 if (rb.
rows () > 0 && rb.
cols () > 0)
620 retval.
xcidx (i) = nz;
629 retval.
xridx (q) = j;
632 assert (
nnz () == retval.
xcidx (nr));
665 return inverse (mattype, info, rcond, 0, 0);
673 return inverse (mattype, info, rcond, 0, 0);
680 return inverse (mattype, info, rcond, 0, 0);
685 double& rcond,
const bool,
686 const bool calccond)
const
694 if (nr == 0 || nc == 0 || nr != nc)
695 (*current_liboctave_error_handler) (
"inverse requires square matrix");
698 int typ = mattype.
type ();
702 (*current_liboctave_error_handler) (
"incorrect matrix type");
718 double tmp = std::abs (v[i]);
735 double& rcond,
const bool,
736 const bool calccond)
const
744 if (nr == 0 || nc == 0 || nr != nc)
745 (*current_liboctave_error_handler) (
"inverse requires square matrix");
748 int typ = mattype.
type ();
753 (*current_liboctave_error_handler) (
"incorrect matrix type");
756 double ainvnorm = 0.;
765 atmp += std::abs (
data (i));
790 retval.
xcidx (i) = cx;
791 retval.
xridx (cx) = i;
792 retval.
xdata (cx) = 1.0;
805 (*current_liboctave_error_handler) (
"division by zero");
810 rpX = retval.
xridx (colXp);
819 v -= retval.
xdata (colXp) *
data (colUp);
824 while (rpX < j && rpU < j && colXp < cx && colUp < nz);
828 colUp =
cidx (j+1) - 1;
832 if (pivot == 0. ||
ridx (colUp) != j)
833 (*current_liboctave_error_handler) (
"division by zero");
843 retval.
xridx (cx) = j;
844 retval.
xdata (cx) = v / pivot;
852 colUp =
cidx (i+1) - 1;
856 if (pivot == 0. ||
ridx (colUp) != i)
857 (*current_liboctave_error_handler) (
"division by zero");
861 retval.
xdata (j) /= pivot;
863 retval.
xcidx (nr) = cx;
908 k <
cidx (jidx+1); k++)
921 (*current_liboctave_error_handler) (
"division by zero");
929 colUp =
cidx (perm[iidx]+1) - 1;
931 colUp =
cidx (perm[iidx]);
935 (*current_liboctave_error_handler) (
"division by zero");
948 nz2 = (2*nz2 < new_cx ? new_cx : 2*nz2);
952 retval.
xcidx (i) = cx;
956 retval.
xridx (cx) = j;
957 retval.
xdata (cx++) = work[j];
961 retval.
xcidx (nr) = cx;
972 i < retval.
cidx (j+1); i++)
973 atmp += std::abs (retval.
data (i));
978 rcond = 1. / ainvnorm / anorm;
986 double& rcond,
bool,
bool calc_cond)
const
990 (*current_liboctave_error_handler)
991 (
"inverse of the null matrix not defined");
994 int typ = mattype.
type (
false);
998 typ = mattype.
type (*
this);
1001 ret = dinverse (mattype, info, rcond,
true, calc_cond);
1003 ret = tinverse (mattype, info, rcond,
true, calc_cond).
transpose ();
1007 ret =
transpose ().tinverse (newtype, info, rcond,
true, calc_cond);
1014 octave::math::sparse_chol<SparseComplexMatrix> fact (*
this, info,
false);
1015 rcond = fact.rcond ();
1021 tinverse (tmp_typ, info, rcond2,
1023 ret =
Q * InvL.
hermitian () * InvL *
Q.transpose ();
1040 octave::math::sparse_lu<SparseComplexMatrix> fact (*
this,
1043 rcond = fact.rcond ();
1051 std::copy_n (
ridx (), nz, ret.
xridx ());
1058 tinverse (tmp_typ, info, rcond2,
1061 tinverse (tmp_typ, info, rcond2,
1063 ret = fact.Pc ().
transpose () * InvU * InvL * fact.Pr ();
1091 #if defined (HAVE_UMFPACK)
1096 if (nr == 0 || nc == 0 || nr != nc)
1105 Matrix Control (UMFPACK_CONTROL, 1);
1109 double tmp = octave::sparse_params::get_key (
"spumoni");
1111 Control (UMFPACK_PRL) = tmp;
1113 tmp = octave::sparse_params::get_key (
"piv_tol");
1116 Control (UMFPACK_SYM_PIVOT_TOLERANCE) = tmp;
1117 Control (UMFPACK_PIVOT_TOLERANCE) = tmp;
1121 tmp = octave::sparse_params::get_key (
"autoamd");
1123 Control (UMFPACK_FIXQ) = tmp;
1126 Control (UMFPACK_SCALE) = UMFPACK_SCALE_NONE;
1137 reinterpret_cast<const double *
> (Ax),
1138 nullptr, 1, control);
1141 Matrix Info (1, UMFPACK_INFO);
1146 reinterpret_cast<const double *
> (Ax),
1147 nullptr,
nullptr, &Symbolic, control, info);
1156 (*current_liboctave_error_handler)
1157 (
"SparseComplexMatrix::determinant symbolic factorization failed");
1167 reinterpret_cast<const double *
> (Ax),
1168 nullptr, Symbolic, &Numeric, control, info);
1171 rcond = Info (UMFPACK_RCOND);
1180 (*current_liboctave_error_handler)
1181 (
"SparseComplexMatrix::determinant numeric factorization failed");
1189 status =
UMFPACK_ZNAME (get_determinant) (c10,
nullptr, &e10,
1197 (*current_liboctave_error_handler)
1198 (
"SparseComplexMatrix::determinant error calculating determinant");
1210 octave_unused_parameter (err);
1211 octave_unused_parameter (rcond);
1213 (*current_liboctave_error_handler)
1214 (
"support for UMFPACK was unavailable or disabled when liboctave was built");
1224 solve_singularity_handler,
bool calc_cond)
const
1233 if (nr != b.
rows ())
1235 (
"matrix dimension mismatch solution of linear equations");
1237 if (nr == 0 || nc == 0 || b.
cols () == 0)
1242 int typ = mattype.
type ();
1246 (*current_liboctave_error_handler) (
"incorrect matrix type");
1252 retval(i, j) = b(i, j) /
data (i);
1257 retval(k, j) = b(
ridx (i), j) /
data (i);
1265 double tmp = std::abs (
data (i));
1271 rcond = dmin / dmax;
1283 solve_singularity_handler,
1284 bool calc_cond)
const
1293 if (nr != b.
rows ())
1295 (
"matrix dimension mismatch solution of linear equations");
1297 if (nr == 0 || nc == 0 || b.
cols () == 0)
1302 int typ = mattype.
type ();
1306 (*current_liboctave_error_handler) (
"incorrect matrix type");
1312 retval.
xcidx (0) = 0;
1319 if (b.
ridx (i) >= nm)
1324 retval.
xcidx (j+1) = ii;
1334 for (k = b.
cidx (j); k < b.
cidx (j+1); k++)
1342 retval.
xridx (ii) = l;
1346 retval.
xcidx (j+1) = ii;
1355 double tmp = std::abs (
data (i));
1361 rcond = dmin / dmax;
1373 solve_singularity_handler,
1374 bool calc_cond)
const
1383 if (nr != b.
rows ())
1385 (
"matrix dimension mismatch solution of linear equations");
1387 if (nr == 0 || nc == 0 || b.
cols () == 0)
1392 int typ = mattype.
type ();
1396 (*current_liboctave_error_handler) (
"incorrect matrix type");
1402 retval(i, j) = b(i, j) /
data (i);
1407 retval(k, j) = b(
ridx (i), j) /
data (i);
1415 double tmp = std::abs (
data (i));
1421 rcond = dmin / dmax;
1433 solve_singularity_handler,
1434 bool calc_cond)
const
1443 if (nr != b.
rows ())
1445 (
"matrix dimension mismatch solution of linear equations");
1447 if (nr == 0 || nc == 0 || b.
cols () == 0)
1452 int typ = mattype.
type ();
1456 (*current_liboctave_error_handler) (
"incorrect matrix type");
1462 retval.
xcidx (0) = 0;
1469 if (b.
ridx (i) >= nm)
1474 retval.
xcidx (j+1) = ii;
1484 for (k = b.
cidx (j); k < b.
cidx (j+1); k++)
1492 retval.
xridx (ii) = l;
1496 retval.
xcidx (j+1) = ii;
1505 double tmp = std::abs (
data (i));
1511 rcond = dmin / dmax;
1523 solve_singularity_handler sing_handler,
1524 bool calc_cond)
const
1533 if (nr != b.
rows ())
1535 (
"matrix dimension mismatch solution of linear equations");
1537 if (nr == 0 || nc == 0 || b.
cols () == 0)
1542 int typ = mattype.
type ();
1546 (*current_liboctave_error_handler) (
"incorrect matrix type");
1549 double ainvnorm = 0.;
1560 atmp += std::abs (
data (i));
1568 retval.
resize (nc, b_nc);
1589 goto triangular_error;
1595 i <
cidx (kidx+1)-1; i++)
1598 work[iidx] = work[iidx] - tmp *
data (i);
1604 retval(perm[i], j) = work[i];
1626 i <
cidx (iidx+1)-1; i++)
1629 work[idx2] = work[idx2] - tmp *
data (i);
1636 atmp += std::abs (work[i]);
1639 if (atmp > ainvnorm)
1642 rcond = 1. / ainvnorm / anorm;
1648 retval.
resize (nc, b_nc);
1665 goto triangular_error;
1673 work[iidx] = work[iidx] - tmp *
data (i);
1679 retval.
xelem (i, j) = work[i];
1699 i <
cidx (k+1)-1; i++)
1702 work[iidx] = work[iidx] - tmp *
data (i);
1709 atmp += std::abs (work[i]);
1712 if (atmp > ainvnorm)
1715 rcond = 1. / ainvnorm / anorm;
1724 sing_handler (rcond);
1731 volatile double rcond_plus_one = rcond + 1.0;
1739 sing_handler (rcond);
1753 solve_singularity_handler sing_handler,
1754 bool calc_cond)
const
1763 if (nr != b.
rows ())
1765 (
"matrix dimension mismatch solution of linear equations");
1767 if (nr == 0 || nc == 0 || b.
cols () == 0)
1772 int typ = mattype.
type ();
1776 (*current_liboctave_error_handler) (
"incorrect matrix type");
1779 double ainvnorm = 0.;
1789 atmp += std::abs (
data (i));
1798 retval.
xcidx (0) = 0;
1828 goto triangular_error;
1834 i <
cidx (kidx+1)-1; i++)
1837 work[iidx] = work[iidx] - tmp *
data (i);
1849 if (ii + new_nnz > x_nz)
1858 if (work[rperm[i]] != 0.)
1860 retval.
xridx (ii) = i;
1861 retval.
xdata (ii++) = work[rperm[i]];
1863 retval.
xcidx (j+1) = ii;
1887 i <
cidx (iidx+1)-1; i++)
1890 work[idx2] = work[idx2] - tmp *
data (i);
1897 atmp += std::abs (work[i]);
1900 if (atmp > ainvnorm)
1903 rcond = 1. / ainvnorm / anorm;
1925 goto triangular_error;
1933 work[iidx] = work[iidx] - tmp *
data (i);
1945 if (ii + new_nnz > x_nz)
1956 retval.
xridx (ii) = i;
1957 retval.
xdata (ii++) = work[i];
1959 retval.
xcidx (j+1) = ii;
1981 i <
cidx (k+1)-1; i++)
1984 work[iidx] = work[iidx] - tmp *
data (i);
1991 atmp += std::abs (work[i]);
1994 if (atmp > ainvnorm)
1997 rcond = 1. / ainvnorm / anorm;
2006 sing_handler (rcond);
2013 volatile double rcond_plus_one = rcond + 1.0;
2021 sing_handler (rcond);
2034 solve_singularity_handler sing_handler,
2035 bool calc_cond)
const
2044 if (nr != b.
rows ())
2046 (
"matrix dimension mismatch solution of linear equations");
2048 if (nr == 0 || nc == 0 || b.
cols () == 0)
2053 int typ = mattype.
type ();
2057 (*current_liboctave_error_handler) (
"incorrect matrix type");
2060 double ainvnorm = 0.;
2071 atmp += std::abs (
data (i));
2079 retval.
resize (nc, b_nc);
2100 goto triangular_error;
2106 i <
cidx (kidx+1)-1; i++)
2109 work[iidx] = work[iidx] - tmp *
data (i);
2115 retval(perm[i], j) = work[i];
2137 i <
cidx (iidx+1)-1; i++)
2140 work[idx2] = work[idx2] - tmp *
data (i);
2147 atmp += std::abs (work[i]);
2150 if (atmp > ainvnorm)
2153 rcond = 1. / ainvnorm / anorm;
2159 retval.
resize (nc, b_nc);
2176 goto triangular_error;
2184 work[iidx] = work[iidx] - tmp *
data (i);
2190 retval.
xelem (i, j) = work[i];
2210 i <
cidx (k+1)-1; i++)
2213 work[iidx] = work[iidx] - tmp *
data (i);
2220 atmp += std::abs (work[i]);
2223 if (atmp > ainvnorm)
2226 rcond = 1. / ainvnorm / anorm;
2235 sing_handler (rcond);
2242 volatile double rcond_plus_one = rcond + 1.0;
2250 sing_handler (rcond);
2264 solve_singularity_handler sing_handler,
2265 bool calc_cond)
const
2274 if (nr != b.
rows ())
2276 (
"matrix dimension mismatch solution of linear equations");
2278 if (nr == 0 || nc == 0 || b.
cols () == 0)
2283 int typ = mattype.
type ();
2287 (*current_liboctave_error_handler) (
"incorrect matrix type");
2290 double ainvnorm = 0.;
2300 atmp += std::abs (
data (i));
2309 retval.
xcidx (0) = 0;
2339 goto triangular_error;
2345 i <
cidx (kidx+1)-1; i++)
2348 work[iidx] = work[iidx] - tmp *
data (i);
2360 if (ii + new_nnz > x_nz)
2369 if (work[rperm[i]] != 0.)
2371 retval.
xridx (ii) = i;
2372 retval.
xdata (ii++) = work[rperm[i]];
2374 retval.
xcidx (j+1) = ii;
2398 i <
cidx (iidx+1)-1; i++)
2401 work[idx2] = work[idx2] - tmp *
data (i);
2408 atmp += std::abs (work[i]);
2411 if (atmp > ainvnorm)
2414 rcond = 1. / ainvnorm / anorm;
2436 goto triangular_error;
2444 work[iidx] = work[iidx] - tmp *
data (i);
2456 if (ii + new_nnz > x_nz)
2467 retval.
xridx (ii) = i;
2468 retval.
xdata (ii++) = work[i];
2470 retval.
xcidx (j+1) = ii;
2492 i <
cidx (k+1)-1; i++)
2495 work[iidx] = work[iidx] - tmp *
data (i);
2502 atmp += std::abs (work[i]);
2505 if (atmp > ainvnorm)
2508 rcond = 1. / ainvnorm / anorm;
2517 sing_handler (rcond);
2524 volatile double rcond_plus_one = rcond + 1.0;
2532 sing_handler (rcond);
2546 solve_singularity_handler sing_handler,
2547 bool calc_cond)
const
2556 if (nr != b.
rows ())
2558 (
"matrix dimension mismatch solution of linear equations");
2560 if (nr == 0 || nc == 0 || b.
cols () == 0)
2565 int typ = mattype.
type ();
2569 (*current_liboctave_error_handler) (
"incorrect matrix type");
2572 double ainvnorm = 0.;
2583 atmp += std::abs (
data (i));
2591 retval.
resize (nc, b_nc);
2600 work[perm[i]] = b(i, j);
2610 if (perm[
ridx (i)] < minr)
2612 minr = perm[
ridx (i)];
2616 if (minr != k ||
data (mini) == 0.)
2619 goto triangular_error;
2630 work[iidx] = work[iidx] - tmp *
data (i);
2636 retval(i, j) = work[i];
2657 i <
cidx (k+1); i++)
2658 if (perm[
ridx (i)] < minr)
2660 minr = perm[
ridx (i)];
2667 i <
cidx (k+1); i++)
2673 work[iidx] = work[iidx] - tmp *
data (i);
2681 atmp += std::abs (work[i]);
2684 if (atmp > ainvnorm)
2687 rcond = 1. / ainvnorm / anorm;
2693 retval.
resize (nc, b_nc, 0.);
2708 goto triangular_error;
2716 work[iidx] = work[iidx] - tmp *
data (i);
2721 retval.
xelem (i, j) = work[i];
2742 i <
cidx (k+1); i++)
2745 work[iidx] = work[iidx] - tmp *
data (i);
2752 atmp += std::abs (work[i]);
2755 if (atmp > ainvnorm)
2758 rcond = 1. / ainvnorm / anorm;
2766 sing_handler (rcond);
2773 volatile double rcond_plus_one = rcond + 1.0;
2781 sing_handler (rcond);
2795 solve_singularity_handler sing_handler,
2796 bool calc_cond)
const
2806 if (nr != b.
rows ())
2808 (
"matrix dimension mismatch solution of linear equations");
2810 if (nr == 0 || nc == 0 || b.
cols () == 0)
2815 int typ = mattype.
type ();
2819 (*current_liboctave_error_handler) (
"incorrect matrix type");
2822 double ainvnorm = 0.;
2832 atmp += std::abs (
data (i));
2841 retval.
xcidx (0) = 0;
2855 work[perm[b.
ridx (i)]] = b.
data (i);
2865 if (perm[
ridx (i)] < minr)
2867 minr = perm[
ridx (i)];
2871 if (minr != k ||
data (mini) == 0.)
2874 goto triangular_error;
2885 work[iidx] = work[iidx] - tmp *
data (i);
2897 if (ii + new_nnz > x_nz)
2908 retval.
xridx (ii) = i;
2909 retval.
xdata (ii++) = work[i];
2911 retval.
xcidx (j+1) = ii;
2934 i <
cidx (k+1); i++)
2935 if (perm[
ridx (i)] < minr)
2937 minr = perm[
ridx (i)];
2944 i <
cidx (k+1); i++)
2950 work[iidx] = work[iidx] - tmp *
data (i);
2958 atmp += std::abs (work[i]);
2961 if (atmp > ainvnorm)
2964 rcond = 1. / ainvnorm / anorm;
2985 goto triangular_error;
2993 work[iidx] = work[iidx] - tmp *
data (i);
3005 if (ii + new_nnz > x_nz)
3016 retval.
xridx (ii) = i;
3017 retval.
xdata (ii++) = work[i];
3019 retval.
xcidx (j+1) = ii;
3042 i <
cidx (k+1); i++)
3045 work[iidx] = work[iidx] - tmp *
data (i);
3052 atmp += std::abs (work[i]);
3055 if (atmp > ainvnorm)
3058 rcond = 1. / ainvnorm / anorm;
3067 sing_handler (rcond);
3074 volatile double rcond_plus_one = rcond + 1.0;
3082 sing_handler (rcond);
3096 solve_singularity_handler sing_handler,
3097 bool calc_cond)
const
3106 if (nr != b.
rows ())
3108 (
"matrix dimension mismatch solution of linear equations");
3110 if (nr == 0 || nc == 0 || b.
cols () == 0)
3115 int typ = mattype.
type ();
3119 (*current_liboctave_error_handler) (
"incorrect matrix type");
3122 double ainvnorm = 0.;
3133 atmp += std::abs (
data (i));
3141 retval.
resize (nc, b_nc);
3150 work[perm[i]] = b(i, j);
3160 if (perm[
ridx (i)] < minr)
3162 minr = perm[
ridx (i)];
3166 if (minr != k ||
data (mini) == 0.)
3169 goto triangular_error;
3180 work[iidx] = work[iidx] - tmp *
data (i);
3186 retval(i, j) = work[i];
3207 i <
cidx (k+1); i++)
3208 if (perm[
ridx (i)] < minr)
3210 minr = perm[
ridx (i)];
3217 i <
cidx (k+1); i++)
3223 work[iidx] = work[iidx] - tmp *
data (i);
3231 atmp += std::abs (work[i]);
3234 if (atmp > ainvnorm)
3237 rcond = 1. / ainvnorm / anorm;
3243 retval.
resize (nc, b_nc, 0.);
3259 goto triangular_error;
3267 work[iidx] = work[iidx] - tmp *
data (i);
3273 retval.
xelem (i, j) = work[i];
3294 i <
cidx (k+1); i++)
3297 work[iidx] = work[iidx] - tmp *
data (i);
3304 atmp += std::abs (work[i]);
3307 if (atmp > ainvnorm)
3310 rcond = 1. / ainvnorm / anorm;
3319 sing_handler (rcond);
3326 volatile double rcond_plus_one = rcond + 1.0;
3334 sing_handler (rcond);
3348 solve_singularity_handler sing_handler,
3349 bool calc_cond)
const
3358 if (nr != b.
rows ())
3360 (
"matrix dimension mismatch solution of linear equations");
3362 if (nr == 0 || nc == 0 || b.
cols () == 0)
3367 int typ = mattype.
type ();
3371 (*current_liboctave_error_handler) (
"incorrect matrix type");
3374 double ainvnorm = 0.;
3384 atmp += std::abs (
data (i));
3393 retval.
xcidx (0) = 0;
3407 work[perm[b.
ridx (i)]] = b.
data (i);
3417 if (perm[
ridx (i)] < minr)
3419 minr = perm[
ridx (i)];
3423 if (minr != k ||
data (mini) == 0.)
3426 goto triangular_error;
3437 work[iidx] = work[iidx] - tmp *
data (i);
3449 if (ii + new_nnz > x_nz)
3460 retval.
xridx (ii) = i;
3461 retval.
xdata (ii++) = work[i];
3463 retval.
xcidx (j+1) = ii;
3486 i <
cidx (k+1); i++)
3487 if (perm[
ridx (i)] < minr)
3489 minr = perm[
ridx (i)];
3496 i <
cidx (k+1); i++)
3502 work[iidx] = work[iidx] - tmp *
data (i);
3510 atmp += std::abs (work[i]);
3513 if (atmp > ainvnorm)
3516 rcond = 1. / ainvnorm / anorm;
3537 goto triangular_error;
3545 work[iidx] = work[iidx] - tmp *
data (i);
3557 if (ii + new_nnz > x_nz)
3568 retval.
xridx (ii) = i;
3569 retval.
xdata (ii++) = work[i];
3571 retval.
xcidx (j+1) = ii;
3594 i <
cidx (k+1); i++)
3597 work[iidx] = work[iidx] - tmp *
data (i);
3604 atmp += std::abs (work[i]);
3607 if (atmp > ainvnorm)
3610 rcond = 1. / ainvnorm / anorm;
3619 sing_handler (rcond);
3626 volatile double rcond_plus_one = rcond + 1.0;
3634 sing_handler (rcond);
3648 solve_singularity_handler sing_handler,
3649 bool calc_cond)
const
3657 if (nr != nc || nr != b.
rows ())
3658 (*current_liboctave_error_handler)
3659 (
"matrix dimension mismatch solution of linear equations");
3661 if (nr == 0 || b.
cols () == 0)
3664 (*current_liboctave_error_handler)
3665 (
"calculation of condition number not implemented");
3669 volatile int typ = mattype.
type ();
3703 else if (
ridx (i) == j + 1)
3714 F77_INT tmp_nr = octave::to_f77_int (nr);
3747 DL[j] =
data (ii++);
3748 DU[j] =
data (ii++);
3750 D[nc-1] =
data (ii);
3767 else if (
ridx (i) == j + 1)
3769 else if (
ridx (i) == j - 1)
3780 F77_INT tmp_nr = octave::to_f77_int (nr);
3797 sing_handler (rcond);
3808 (*current_liboctave_error_handler) (
"incorrect matrix type");
3817 solve_singularity_handler sing_handler,
3818 bool calc_cond)
const
3826 if (nr != nc || nr != b.
rows ())
3827 (*current_liboctave_error_handler)
3828 (
"matrix dimension mismatch solution of linear equations");
3830 if (nr == 0 || b.
cols () == 0)
3833 (*current_liboctave_error_handler)
3834 (
"calculation of condition number not implemented");
3838 int typ = mattype.
type ();
3850 F77_INT *pipvt = ipvt.fortran_vec ();
3859 DL[j] =
data (ii++);
3860 DU[j] =
data (ii++);
3862 D[nc-1] =
data (ii);
3879 else if (
ridx (i) == j + 1)
3881 else if (
ridx (i) == j - 1)
3886 F77_INT tmp_nr = octave::to_f77_int (nr);
3905 sing_handler (rcond);
3918 retval.
xcidx (0) = 0;
3932 (F77_CONST_CHAR_ARG2 (&job, 1),
3938 F77_CHAR_ARG_LEN (1)));
3949 if (ii + new_nnz > x_nz)
3960 retval.
xridx (ii) = i;
3961 retval.
xdata (ii++) = work[i];
3963 retval.
xcidx (j+1) = ii;
3970 (*current_liboctave_error_handler) (
"incorrect matrix type");
3979 solve_singularity_handler sing_handler,
3980 bool calc_cond)
const
3988 if (nr != nc || nr != b.
rows ())
3989 (*current_liboctave_error_handler)
3990 (
"matrix dimension mismatch solution of linear equations");
3992 if (nr == 0 || b.
cols () == 0)
3995 (*current_liboctave_error_handler)
3996 (
"calculation of condition number not implemented");
4000 volatile int typ = mattype.
type ();
4034 else if (
ridx (i) == j + 1)
4047 F77_INT tmp_nr = octave::to_f77_int (nr);
4078 DL[j] =
data (ii++);
4079 DU[j] =
data (ii++);
4081 D[nc-1] =
data (ii);
4098 else if (
ridx (i) == j + 1)
4100 else if (
ridx (i) == j - 1)
4113 F77_INT tmp_nr = octave::to_f77_int (nr);
4130 sing_handler (rcond);
4138 (*current_liboctave_error_handler) (
"incorrect matrix type");
4145 SparseComplexMatrix::trisolve (
MatrixType& mattype,
4148 solve_singularity_handler sing_handler,
4149 bool calc_cond)
const
4157 if (nr != nc || nr != b.
rows ())
4158 (*current_liboctave_error_handler)
4159 (
"matrix dimension mismatch solution of linear equations");
4161 if (nr == 0 || b.
cols () == 0)
4164 (*current_liboctave_error_handler)
4165 (
"calculation of condition number not implemented");
4169 int typ = mattype.
type ();
4181 F77_INT *pipvt = ipvt.fortran_vec ();
4190 DL[j] =
data (ii++);
4191 DU[j] =
data (ii++);
4193 D[nc-1] =
data (ii);
4210 else if (
ridx (i) == j + 1)
4212 else if (
ridx (i) == j - 1)
4217 F77_INT tmp_nr = octave::to_f77_int (nr);
4236 sing_handler (rcond);
4256 retval.
xcidx (0) = 0;
4260 for (
F77_INT i = 0; i < b_nr; i++)
4264 (F77_CONST_CHAR_ARG2 (&job, 1),
4270 F77_CHAR_ARG_LEN (1)));
4278 (*current_liboctave_error_handler)
4279 (
"SparseComplexMatrix::solve solve failed");
4292 if (ii + new_nnz > x_nz)
4303 retval.
xridx (ii) = i;
4304 retval.
xdata (ii++) = Bx[i];
4307 retval.
xcidx (j+1) = ii;
4314 (*current_liboctave_error_handler) (
"incorrect matrix type");
4323 solve_singularity_handler sing_handler,
4324 bool calc_cond)
const
4332 if (nr != nc || nr != b.
rows ())
4333 (*current_liboctave_error_handler)
4334 (
"matrix dimension mismatch solution of linear equations");
4336 if (nr == 0 || b.
cols () == 0)
4341 volatile int typ = mattype.
type ();
4349 Complex *tmp_data = m_band.fortran_vec ();
4355 for (
F77_INT j = 0; j < ldm; j++)
4357 tmp_data[ii++] = 0.;
4365 m_band(ri - j, j) =
data (i);
4371 anorm = m_band.abs ().sum ().row (0).max ();
4373 F77_INT tmp_nr = octave::to_f77_int (nr);
4378 F77_XFCN (zpbtrf, ZPBTRF, (F77_CONST_CHAR_ARG2 (&job, 1),
4380 F77_CHAR_ARG_LEN (1)));
4398 Complex *pz = z.fortran_vec ();
4400 double *piz = iz.fortran_vec ();
4403 (F77_CONST_CHAR_ARG2 (&job, 1),
4406 F77_CHAR_ARG_LEN (1)));
4413 volatile double rcond_plus_one = rcond + 1.0;
4421 sing_handler (rcond);
4440 (F77_CONST_CHAR_ARG2 (&job, 1),
4443 F77_CHAR_ARG_LEN (1)));
4450 (*current_liboctave_error_handler)
4451 (
"SparseMatrix::solve solve failed");
4463 F77_INT ldm = n_upper + 2 * n_lower + 1;
4466 Complex *tmp_data = m_band.fortran_vec ();
4472 for (
F77_INT j = 0; j < ldm; j++)
4474 tmp_data[ii++] = 0.;
4479 m_band(
ridx (i) - j + n_lower + n_upper, j) =
data (i);
4489 atmp += std::abs (
data (i));
4496 F77_INT *pipvt = ipvt.fortran_vec ();
4498 F77_INT tmp_nr = octave::to_f77_int (nr);
4502 F77_XFCN (zgbtrf, ZGBTRF, (tmp_nr, tmp_nr, n_lower, n_upper,
4504 ldm, pipvt, tmp_err));
4517 sing_handler (rcond);
4529 Complex *pz = z.fortran_vec ();
4531 double *piz = iz.fortran_vec ();
4533 F77_INT tmp_nc = octave::to_f77_int (nc);
4536 (F77_CONST_CHAR_ARG2 (&job, 1),
4539 F77_CHAR_ARG_LEN (1)));
4546 volatile double rcond_plus_one = rcond + 1.0;
4554 sing_handler (rcond);
4574 (F77_CONST_CHAR_ARG2 (&job, 1),
4577 F77_CHAR_ARG_LEN (1)));
4584 (*current_liboctave_error_handler) (
"incorrect matrix type");
4593 solve_singularity_handler sing_handler,
4594 bool calc_cond)
const
4602 if (nr != nc || nr != b.
rows ())
4603 (*current_liboctave_error_handler)
4604 (
"matrix dimension mismatch solution of linear equations");
4606 if (nr == 0 || b.
cols () == 0)
4611 volatile int typ = mattype.
type ();
4620 Complex *tmp_data = m_band.fortran_vec ();
4626 for (
F77_INT j = 0; j < ldm; j++)
4628 tmp_data[ii++] = 0.;
4636 m_band(ri - j, j) =
data (i);
4642 anorm = m_band.abs ().sum ().row (0).max ();
4644 F77_INT tmp_nr = octave::to_f77_int (nr);
4649 F77_XFCN (zpbtrf, ZPBTRF, (F77_CONST_CHAR_ARG2 (&job, 1),
4651 F77_CHAR_ARG_LEN (1)));
4667 Complex *pz = z.fortran_vec ();
4669 double *piz = iz.fortran_vec ();
4672 (F77_CONST_CHAR_ARG2 (&job, 1),
4675 F77_CHAR_ARG_LEN (1)));
4682 volatile double rcond_plus_one = rcond + 1.0;
4690 sing_handler (rcond);
4712 retval.
xcidx (0) = 0;
4715 for (
F77_INT i = 0; i < b_nr; i++)
4716 Bx[i] = b.
elem (i, j);
4719 (F77_CONST_CHAR_ARG2 (&job, 1),
4722 F77_CHAR_ARG_LEN (1)));
4729 (*current_liboctave_error_handler)
4730 (
"SparseComplexMatrix::solve solve failed");
4744 sz = (
static_cast<double> (b_nc) - j) / b_nc
4746 sz = x_nz + (sz > 100 ? sz : 100);
4750 retval.
xdata (ii) = tmp;
4751 retval.
xridx (ii++) = i;
4754 retval.
xcidx (j+1) = ii;
4767 F77_INT ldm = n_upper + 2 * n_lower + 1;
4770 Complex *tmp_data = m_band.fortran_vec ();
4776 for (
F77_INT j = 0; j < ldm; j++)
4778 tmp_data[ii++] = 0.;
4783 m_band(
ridx (i) - j + n_lower + n_upper, j) =
data (i);
4793 atmp += std::abs (
data (i));
4800 F77_INT *pipvt = ipvt.fortran_vec ();
4802 F77_INT tmp_nr = octave::to_f77_int (nr);
4806 F77_XFCN (zgbtrf, ZGBTRF, (tmp_nr, tmp_nr, n_lower, n_upper,
4808 ldm, pipvt, tmp_err));
4819 sing_handler (rcond);
4831 Complex *pz = z.fortran_vec ();
4833 double *piz = iz.fortran_vec ();
4835 F77_INT tmp_nc = octave::to_f77_int (nc);
4838 (F77_CONST_CHAR_ARG2 (&job, 1),
4841 F77_CHAR_ARG_LEN (1)));
4848 volatile double rcond_plus_one = rcond + 1.0;
4856 sing_handler (rcond);
4873 retval.
xcidx (0) = 0;
4883 i < b.
cidx (j+1); i++)
4887 (F77_CONST_CHAR_ARG2 (&job, 1),
4890 F77_CHAR_ARG_LEN (1)));
4901 if (ii + new_nnz > x_nz)
4912 retval.
xridx (ii) = i;
4913 retval.
xdata (ii++) = work[i];
4915 retval.
xcidx (j+1) = ii;
4923 (*current_liboctave_error_handler) (
"incorrect matrix type");
4932 solve_singularity_handler sing_handler,
4933 bool calc_cond)
const
4941 if (nr != nc || nr != b.
rows ())
4942 (*current_liboctave_error_handler)
4943 (
"matrix dimension mismatch solution of linear equations");
4945 if (nr == 0 || b.
cols () == 0)
4950 volatile int typ = mattype.
type ();
4959 Complex *tmp_data = m_band.fortran_vec ();
4965 for (
F77_INT j = 0; j < ldm; j++)
4967 tmp_data[ii++] = 0.;
4975 m_band(ri - j, j) =
data (i);
4981 anorm = m_band.abs ().sum ().row (0).max ();
4983 F77_INT tmp_nr = octave::to_f77_int (nr);
4988 F77_XFCN (zpbtrf, ZPBTRF, (F77_CONST_CHAR_ARG2 (&job, 1),
4990 F77_CHAR_ARG_LEN (1)));
5008 Complex *pz = z.fortran_vec ();
5010 double *piz = iz.fortran_vec ();
5013 (F77_CONST_CHAR_ARG2 (&job, 1),
5016 F77_CHAR_ARG_LEN (1)));
5023 volatile double rcond_plus_one = rcond + 1.0;
5031 sing_handler (rcond);
5049 (F77_CONST_CHAR_ARG2 (&job, 1),
5052 F77_CHAR_ARG_LEN (1)));
5059 (*current_liboctave_error_handler)
5060 (
"SparseComplexMatrix::solve solve failed");
5072 F77_INT ldm = n_upper + 2 * n_lower + 1;
5075 Complex *tmp_data = m_band.fortran_vec ();
5081 for (
F77_INT j = 0; j < ldm; j++)
5083 tmp_data[ii++] = 0.;
5088 m_band(
ridx (i) - j + n_lower + n_upper, j) =
data (i);
5098 atmp += std::abs (
data (i));
5105 F77_INT *pipvt = ipvt.fortran_vec ();
5107 F77_INT tmp_nr = octave::to_f77_int (nr);
5111 F77_XFCN (zgbtrf, ZGBTRF, (tmp_nr, tmp_nr, n_lower, n_upper,
5113 ldm, pipvt, tmp_err));
5124 sing_handler (rcond);
5136 Complex *pz = z.fortran_vec ();
5138 double *piz = iz.fortran_vec ();
5140 F77_INT tmp_nc = octave::to_f77_int (nc);
5143 (F77_CONST_CHAR_ARG2 (&job, 1),
5146 F77_CHAR_ARG_LEN (1)));
5153 volatile double rcond_plus_one = rcond + 1.0;
5161 sing_handler (rcond);
5180 (F77_CONST_CHAR_ARG2 (&job, 1),
5183 F77_CHAR_ARG_LEN (1)));
5190 (*current_liboctave_error_handler) (
"incorrect matrix type");
5199 solve_singularity_handler sing_handler,
5200 bool calc_cond)
const
5208 if (nr != nc || nr != b.
rows ())
5209 (*current_liboctave_error_handler)
5210 (
"matrix dimension mismatch solution of linear equations");
5212 if (nr == 0 || b.
cols () == 0)
5217 volatile int typ = mattype.
type ();
5226 Complex *tmp_data = m_band.fortran_vec ();
5232 for (
F77_INT j = 0; j < ldm; j++)
5234 tmp_data[ii++] = 0.;
5242 m_band(ri - j, j) =
data (i);
5248 anorm = m_band.abs ().sum ().row (0).max ();
5250 F77_INT tmp_nr = octave::to_f77_int (nr);
5255 F77_XFCN (zpbtrf, ZPBTRF, (F77_CONST_CHAR_ARG2 (&job, 1),
5257 F77_CHAR_ARG_LEN (1)));
5276 Complex *pz = z.fortran_vec ();
5278 double *piz = iz.fortran_vec ();
5281 (F77_CONST_CHAR_ARG2 (&job, 1),
5284 F77_CHAR_ARG_LEN (1)));
5291 volatile double rcond_plus_one = rcond + 1.0;
5299 sing_handler (rcond);
5321 retval.
xcidx (0) = 0;
5325 for (
F77_INT i = 0; i < b_nr; i++)
5329 (F77_CONST_CHAR_ARG2 (&job, 1),
5332 F77_CHAR_ARG_LEN (1)));
5339 (*current_liboctave_error_handler)
5340 (
"SparseMatrix::solve solve failed");
5352 if (ii + new_nnz > x_nz)
5363 retval.
xridx (ii) = i;
5364 retval.
xdata (ii++) = Bx[i];
5367 retval.
xcidx (j+1) = ii;
5380 F77_INT ldm = n_upper + 2 * n_lower + 1;
5383 Complex *tmp_data = m_band.fortran_vec ();
5389 for (
F77_INT j = 0; j < ldm; j++)
5391 tmp_data[ii++] = 0.;
5396 m_band(
ridx (i) - j + n_lower + n_upper, j) =
data (i);
5406 atmp += std::abs (
data (i));
5413 F77_INT *pipvt = ipvt.fortran_vec ();
5415 F77_INT tmp_nr = octave::to_f77_int (nr);
5419 F77_XFCN (zgbtrf, ZGBTRF, (tmp_nr, tmp_nr, n_lower, n_upper,
5421 ldm, pipvt, tmp_err));
5432 sing_handler (rcond);
5444 Complex *pz = z.fortran_vec ();
5446 double *piz = iz.fortran_vec ();
5448 F77_INT tmp_nc = octave::to_f77_int (nc);
5451 (F77_CONST_CHAR_ARG2 (&job, 1),
5454 F77_CHAR_ARG_LEN (1)));
5461 volatile double rcond_plus_one = rcond + 1.0;
5469 sing_handler (rcond);
5486 retval.
xcidx (0) = 0;
5497 i < b.
cidx (j+1); i++)
5501 (F77_CONST_CHAR_ARG2 (&job, 1),
5504 F77_CHAR_ARG_LEN (1)));
5515 if (ii + new_nnz > x_nz)
5526 retval.
xridx (ii) = i;
5527 retval.
xdata (ii++) = Bx[i];
5529 retval.
xcidx (j+1) = ii;
5537 (*current_liboctave_error_handler) (
"incorrect matrix type");
5546 solve_singularity_handler sing_handler,
5547 bool calc_cond)
const
5550 void *Numeric =
nullptr;
5553 #if defined (HAVE_UMFPACK)
5556 Control =
Matrix (UMFPACK_CONTROL, 1);
5560 double tmp = octave::sparse_params::get_key (
"spumoni");
5562 Control (UMFPACK_PRL) = tmp;
5563 tmp = octave::sparse_params::get_key (
"piv_tol");
5566 Control (UMFPACK_SYM_PIVOT_TOLERANCE) = tmp;
5567 Control (UMFPACK_PIVOT_TOLERANCE) = tmp;
5571 tmp = octave::sparse_params::get_key (
"autoamd");
5573 Control (UMFPACK_FIXQ) = tmp;
5586 reinterpret_cast<const double *
> (Ax),
5587 nullptr, 1, control);
5590 Info =
Matrix (1, UMFPACK_INFO);
5595 reinterpret_cast<const double *
> (Ax),
5596 nullptr,
nullptr, &Symbolic, control, info);
5606 (*current_liboctave_error_handler)
5607 (
"SparseComplexMatrix::solve symbolic factorization failed");
5616 reinterpret_cast<const double *
> (Ax),
5617 nullptr, Symbolic, &Numeric, control, info);
5621 rcond = Info (UMFPACK_RCOND);
5624 volatile double rcond_plus_one = rcond + 1.0;
5626 if (status == UMFPACK_WARNING_singular_matrix
5634 sing_handler (rcond);
5638 else if (status < 0)
5644 (*current_liboctave_error_handler)
5645 (
"SparseComplexMatrix::solve numeric factorization failed");
5660 octave_unused_parameter (rcond);
5661 octave_unused_parameter (Control);
5662 octave_unused_parameter (Info);
5663 octave_unused_parameter (sing_handler);
5664 octave_unused_parameter (calc_cond);
5666 (*current_liboctave_error_handler)
5667 (
"support for UMFPACK was unavailable or disabled when liboctave was built");
5677 solve_singularity_handler sing_handler,
5678 bool calc_cond)
const
5686 if (nr != nc || nr != b.
rows ())
5687 (*current_liboctave_error_handler)
5688 (
"matrix dimension mismatch solution of linear equations");
5690 if (nr == 0 || b.
cols () == 0)
5695 volatile int typ = mattype.
type ();
5700 #if defined (HAVE_CHOLMOD)
5701 cholmod_common Common;
5702 cholmod_common *cm = &Common;
5706 cm->prefer_zomplex =
false;
5708 double spu = octave::sparse_params::get_key (
"spumoni");
5712 SUITESPARSE_ASSIGN_FPTR (printf_func, cm->print_function,
nullptr);
5716 cm->print =
static_cast<int> (spu) + 2;
5717 SUITESPARSE_ASSIGN_FPTR (printf_func, cm->print_function, &
SparseCholPrint);
5721 SUITESPARSE_ASSIGN_FPTR2 (divcomplex_func, cm->complex_divide, divcomplex);
5722 SUITESPARSE_ASSIGN_FPTR2 (hypot_func, cm->hypotenuse, hypot);
5724 cm->final_ll =
true;
5726 cholmod_sparse Astore;
5727 cholmod_sparse *
A = &Astore;
5737 #if defined (OCTAVE_ENABLE_64)
5738 A->itype = CHOLMOD_LONG;
5740 A->itype = CHOLMOD_INT;
5742 A->dtype = CHOLMOD_DOUBLE;
5744 A->xtype = CHOLMOD_COMPLEX;
5748 cholmod_dense Bstore;
5749 cholmod_dense *
B = &Bstore;
5750 B->nrow = b.
rows ();
5751 B->ncol = b.
cols ();
5753 B->nzmax =
B->nrow *
B->ncol;
5754 B->dtype = CHOLMOD_DOUBLE;
5755 B->xtype = CHOLMOD_REAL;
5757 B->x =
const_cast<double *
> (b.
data ());
5774 volatile double rcond_plus_one = rcond + 1.0;
5782 sing_handler (rcond);
5798 retval.
xelem (i, j) =
static_cast<Complex *
> (X->x)[jr + i];
5804 static char blank_name[] =
" ";
5808 (*current_liboctave_warning_with_id_handler)
5809 (
"Octave:missing-dependency",
5810 "support for CHOLMOD was unavailable or disabled "
5811 "when liboctave was built");
5820 #if defined (HAVE_UMFPACK)
5822 void *Numeric = factorize (err, rcond, Control, Info,
5823 sing_handler, calc_cond);
5828 Control (UMFPACK_IRSTEP) = 1;
5837 #if defined (UMFPACK_SEPARATE_SPLIT)
5838 const double *Bx = b.
data ();
5845 retval.
resize (b_nr, b_nc);
5850 #if defined (UMFPACK_SEPARATE_SPLIT)
5854 reinterpret_cast<const double *
> (Ax),
5856 reinterpret_cast<double *
> (&Xx[iidx]),
5858 &Bx[iidx], Bz, Numeric,
5862 Bz[i] = b.
elem (i, j);
5867 reinterpret_cast<const double *
> (Ax),
5869 reinterpret_cast<double *
> (&Xx[iidx]),
5871 reinterpret_cast<const double *
> (Bz),
5881 (*current_liboctave_error_handler)
5882 (
"SparseComplexMatrix::solve solve failed");
5897 octave_unused_parameter (rcond);
5898 octave_unused_parameter (sing_handler);
5899 octave_unused_parameter (calc_cond);
5901 (*current_liboctave_error_handler)
5902 (
"support for UMFPACK was unavailable or disabled "
5903 "when liboctave was built");
5907 (*current_liboctave_error_handler) (
"incorrect matrix type");
5916 solve_singularity_handler sing_handler,
5917 bool calc_cond)
const
5925 if (nr != nc || nr != b.
rows ())
5926 (*current_liboctave_error_handler)
5927 (
"matrix dimension mismatch solution of linear equations");
5929 if (nr == 0 || b.
cols () == 0)
5934 volatile int typ = mattype.
type ();
5939 #if defined (HAVE_CHOLMOD)
5940 cholmod_common Common;
5941 cholmod_common *cm = &Common;
5945 cm->prefer_zomplex =
false;
5947 double spu = octave::sparse_params::get_key (
"spumoni");
5951 SUITESPARSE_ASSIGN_FPTR (printf_func, cm->print_function,
nullptr);
5955 cm->print =
static_cast<int> (spu) + 2;
5956 SUITESPARSE_ASSIGN_FPTR (printf_func, cm->print_function, &
SparseCholPrint);
5960 SUITESPARSE_ASSIGN_FPTR2 (divcomplex_func, cm->complex_divide, divcomplex);
5961 SUITESPARSE_ASSIGN_FPTR2 (hypot_func, cm->hypotenuse, hypot);
5963 cm->final_ll =
true;
5965 cholmod_sparse Astore;
5966 cholmod_sparse *
A = &Astore;
5976 #if defined (OCTAVE_ENABLE_64)
5977 A->itype = CHOLMOD_LONG;
5979 A->itype = CHOLMOD_INT;
5981 A->dtype = CHOLMOD_DOUBLE;
5983 A->xtype = CHOLMOD_COMPLEX;
5987 cholmod_sparse Bstore;
5988 cholmod_sparse *
B = &Bstore;
5989 B->nrow = b.
rows ();
5990 B->ncol = b.
cols ();
5993 B->nzmax = b.
nnz ();
5997 #if defined (OCTAVE_ENABLE_64)
5998 B->itype = CHOLMOD_LONG;
6000 B->itype = CHOLMOD_INT;
6002 B->dtype = CHOLMOD_DOUBLE;
6004 B->xtype = CHOLMOD_REAL;
6023 volatile double rcond_plus_one = rcond + 1.0;
6031 sing_handler (rcond);
6040 cholmod_sparse *X =
CHOLMOD_NAME(spsolve) (CHOLMOD_A, L,
B, cm);
6048 j <= static_cast<octave_idx_type> (X->ncol); j++)
6061 static char blank_name[] =
" ";
6065 (*current_liboctave_warning_with_id_handler)
6066 (
"Octave:missing-dependency",
6067 "support for CHOLMOD was unavailable or disabled "
6068 "when liboctave was built");
6077 #if defined (HAVE_UMFPACK)
6079 void *Numeric = factorize (err, rcond, Control, Info,
6080 sing_handler, calc_cond);
6085 Control (UMFPACK_IRSTEP) = 1;
6095 #if defined (UMFPACK_SEPARATE_SPLIT)
6112 retval.
xcidx (0) = 0;
6116 #if defined (UMFPACK_SEPARATE_SPLIT)
6118 Bx[i] = b.
elem (i, j);
6123 reinterpret_cast<const double *
> (Ax),
6125 reinterpret_cast<double *
> (Xx),
6127 Bx, Bz, Numeric, control,
6131 Bz[i] = b.
elem (i, j);
6136 reinterpret_cast<const double *
> (Ax),
6138 reinterpret_cast<double *
> (Xx),
6140 reinterpret_cast<double *
> (Bz),
6150 (*current_liboctave_error_handler)
6151 (
"SparseComplexMatrix::solve solve failed");
6166 sz = (
static_cast<double> (b_nc) - j) / b_nc
6168 sz = x_nz + (sz > 100 ? sz : 100);
6172 retval.
xdata (ii) = tmp;
6173 retval.
xridx (ii++) = i;
6176 retval.
xcidx (j+1) = ii;
6189 octave_unused_parameter (rcond);
6190 octave_unused_parameter (sing_handler);
6191 octave_unused_parameter (calc_cond);
6193 (*current_liboctave_error_handler)
6194 (
"support for UMFPACK was unavailable or disabled "
6195 "when liboctave was built");
6199 (*current_liboctave_error_handler) (
"incorrect matrix type");
6208 solve_singularity_handler sing_handler,
6209 bool calc_cond)
const
6217 if (nr != nc || nr != b.
rows ())
6218 (*current_liboctave_error_handler)
6219 (
"matrix dimension mismatch solution of linear equations");
6221 if (nr == 0 || b.
cols () == 0)
6226 volatile int typ = mattype.
type ();
6231 #if defined (HAVE_CHOLMOD)
6232 cholmod_common Common;
6233 cholmod_common *cm = &Common;
6237 cm->prefer_zomplex =
false;
6239 double spu = octave::sparse_params::get_key (
"spumoni");
6243 SUITESPARSE_ASSIGN_FPTR (printf_func, cm->print_function,
nullptr);
6247 cm->print =
static_cast<int> (spu) + 2;
6248 SUITESPARSE_ASSIGN_FPTR (printf_func, cm->print_function, &
SparseCholPrint);
6252 SUITESPARSE_ASSIGN_FPTR2 (divcomplex_func, cm->complex_divide, divcomplex);
6253 SUITESPARSE_ASSIGN_FPTR2 (hypot_func, cm->hypotenuse, hypot);
6255 cm->final_ll =
true;
6257 cholmod_sparse Astore;
6258 cholmod_sparse *
A = &Astore;
6268 #if defined (OCTAVE_ENABLE_64)
6269 A->itype = CHOLMOD_LONG;
6271 A->itype = CHOLMOD_INT;
6273 A->dtype = CHOLMOD_DOUBLE;
6275 A->xtype = CHOLMOD_COMPLEX;
6279 cholmod_dense Bstore;
6280 cholmod_dense *
B = &Bstore;
6281 B->nrow = b.
rows ();
6282 B->ncol = b.
cols ();
6284 B->nzmax =
B->nrow *
B->ncol;
6285 B->dtype = CHOLMOD_DOUBLE;
6286 B->xtype = CHOLMOD_COMPLEX;
6305 volatile double rcond_plus_one = rcond + 1.0;
6313 sing_handler (rcond);
6329 retval.
xelem (i, j) =
static_cast<Complex *
> (X->x)[jr + i];
6335 static char blank_name[] =
" ";
6339 (*current_liboctave_warning_with_id_handler)
6340 (
"Octave:missing-dependency",
6341 "support for CHOLMOD was unavailable or disabled "
6342 "when liboctave was built");
6351 #if defined (HAVE_UMFPACK)
6353 void *Numeric = factorize (err, rcond, Control, Info,
6354 sing_handler, calc_cond);
6359 Control (UMFPACK_IRSTEP) = 1;
6370 retval.
resize (b_nr, b_nc);
6379 reinterpret_cast<const double *
> (Ax),
6381 reinterpret_cast<double *
> (&Xx[iidx]),
6383 reinterpret_cast<const double *
> (&Bx[iidx]),
6384 nullptr, Numeric, control, info);
6391 (*current_liboctave_error_handler)
6392 (
"SparseComplexMatrix::solve solve failed");
6407 octave_unused_parameter (rcond);
6408 octave_unused_parameter (sing_handler);
6409 octave_unused_parameter (calc_cond);
6411 (*current_liboctave_error_handler)
6412 (
"support for UMFPACK was unavailable or disabled "
6413 "when liboctave was built");
6417 (*current_liboctave_error_handler) (
"incorrect matrix type");
6426 solve_singularity_handler sing_handler,
6427 bool calc_cond)
const
6435 if (nr != nc || nr != b.
rows ())
6436 (*current_liboctave_error_handler)
6437 (
"matrix dimension mismatch solution of linear equations");
6439 if (nr == 0 || b.
cols () == 0)
6444 volatile int typ = mattype.
type ();
6449 #if defined (HAVE_CHOLMOD)
6450 cholmod_common Common;
6451 cholmod_common *cm = &Common;
6455 cm->prefer_zomplex =
false;
6457 double spu = octave::sparse_params::get_key (
"spumoni");
6461 SUITESPARSE_ASSIGN_FPTR (printf_func, cm->print_function,
nullptr);
6465 cm->print =
static_cast<int> (spu) + 2;
6466 SUITESPARSE_ASSIGN_FPTR (printf_func, cm->print_function, &
SparseCholPrint);
6470 SUITESPARSE_ASSIGN_FPTR2 (divcomplex_func, cm->complex_divide, divcomplex);
6471 SUITESPARSE_ASSIGN_FPTR2 (hypot_func, cm->hypotenuse, hypot);
6473 cm->final_ll =
true;
6475 cholmod_sparse Astore;
6476 cholmod_sparse *
A = &Astore;
6486 #if defined (OCTAVE_ENABLE_64)
6487 A->itype = CHOLMOD_LONG;
6489 A->itype = CHOLMOD_INT;
6491 A->dtype = CHOLMOD_DOUBLE;
6493 A->xtype = CHOLMOD_COMPLEX;
6497 cholmod_sparse Bstore;
6498 cholmod_sparse *
B = &Bstore;
6499 B->nrow = b.
rows ();
6500 B->ncol = b.
cols ();
6503 B->nzmax = b.
nnz ();
6507 #if defined (OCTAVE_ENABLE_64)
6508 B->itype = CHOLMOD_LONG;
6510 B->itype = CHOLMOD_INT;
6512 B->dtype = CHOLMOD_DOUBLE;
6514 B->xtype = CHOLMOD_COMPLEX;
6533 volatile double rcond_plus_one = rcond + 1.0;
6541 sing_handler (rcond);
6550 cholmod_sparse *X =
CHOLMOD_NAME(spsolve) (CHOLMOD_A, L,
B, cm);
6558 j <= static_cast<octave_idx_type> (X->ncol); j++)
6571 static char blank_name[] =
" ";
6575 (*current_liboctave_warning_with_id_handler)
6576 (
"Octave:missing-dependency",
6577 "support for CHOLMOD was unavailable or disabled "
6578 "when liboctave was built");
6587 #if defined (HAVE_UMFPACK)
6589 void *Numeric = factorize (err, rcond, Control, Info,
6590 sing_handler, calc_cond);
6595 Control (UMFPACK_IRSTEP) = 1;
6615 retval.
xcidx (0) = 0;
6624 reinterpret_cast<const double *
> (Ax),
6626 reinterpret_cast<double *
> (Xx),
6628 reinterpret_cast<double *
> (Bx),
6629 nullptr, Numeric, control, info);
6636 (*current_liboctave_error_handler)
6637 (
"SparseComplexMatrix::solve solve failed");
6652 sz = (
static_cast<double> (b_nc) - j) / b_nc
6654 sz = x_nz + (sz > 100 ? sz : 100);
6658 retval.
xdata (ii) = tmp;
6659 retval.
xridx (ii++) = i;
6662 retval.
xcidx (j+1) = ii;
6667 rcond = Info (UMFPACK_RCOND);
6668 volatile double rcond_plus_one = rcond + 1.0;
6670 if (status == UMFPACK_WARNING_singular_matrix
6676 sing_handler (rcond);
6689 octave_unused_parameter (rcond);
6690 octave_unused_parameter (sing_handler);
6691 octave_unused_parameter (calc_cond);
6693 (*current_liboctave_error_handler)
6694 (
"support for UMFPACK was unavailable or disabled "
6695 "when liboctave was built");
6699 (*current_liboctave_error_handler) (
"incorrect matrix type");
6710 return solve (mattype, b, info, rcond,
nullptr);
6718 return solve (mattype, b, info, rcond,
nullptr);
6725 return solve (mattype, b, info, rcond,
nullptr);
6731 solve_singularity_handler sing_handler,
6732 bool singular_fallback)
const
6735 int typ = mattype.
type (
false);
6738 typ = mattype.
type (*
this);
6741 retval = dsolve (mattype, b, err, rcond, sing_handler,
false);
6743 retval = utsolve (mattype, b, err, rcond, sing_handler,
false);
6745 retval = ltsolve (mattype, b, err, rcond, sing_handler,
false);
6747 retval = bsolve (mattype, b, err, rcond, sing_handler,
false);
6750 retval = trisolve (mattype, b, err, rcond, sing_handler,
false);
6752 retval = fsolve (mattype, b, err, rcond, sing_handler,
true);
6754 (*current_liboctave_error_handler) (
"unknown matrix type");
6759 #if defined (USE_QRSOLVE)
6760 retval =
qrsolve (*
this, b, err);
6775 return solve (mattype, b, info, rcond,
nullptr);
6783 return solve (mattype, b, info, rcond,
nullptr);
6790 return solve (mattype, b, info, rcond,
nullptr);
6796 solve_singularity_handler sing_handler,
6797 bool singular_fallback)
const
6800 int typ = mattype.
type (
false);
6803 typ = mattype.
type (*
this);
6806 retval = dsolve (mattype, b, err, rcond, sing_handler,
false);
6808 retval = utsolve (mattype, b, err, rcond, sing_handler,
false);
6810 retval = ltsolve (mattype, b, err, rcond, sing_handler,
false);
6812 retval = bsolve (mattype, b, err, rcond, sing_handler,
false);
6815 retval = trisolve (mattype, b, err, rcond, sing_handler,
false);
6817 retval = fsolve (mattype, b, err, rcond, sing_handler,
true);
6819 (*current_liboctave_error_handler) (
"unknown matrix type");
6824 #if defined (USE_QRSOLVE)
6825 retval =
qrsolve (*
this, b, err);
6840 return solve (mattype, b, info, rcond,
nullptr);
6848 return solve (mattype, b, info, rcond,
nullptr);
6855 return solve (mattype, b, info, rcond,
nullptr);
6861 solve_singularity_handler sing_handler,
6862 bool singular_fallback)
const
6865 int typ = mattype.
type (
false);
6868 typ = mattype.
type (*
this);
6871 retval = dsolve (mattype, b, err, rcond, sing_handler,
false);
6873 retval = utsolve (mattype, b, err, rcond, sing_handler,
false);
6875 retval = ltsolve (mattype, b, err, rcond, sing_handler,
false);
6877 retval = bsolve (mattype, b, err, rcond, sing_handler,
false);
6880 retval = trisolve (mattype, b, err, rcond, sing_handler,
false);
6882 retval = fsolve (mattype, b, err, rcond, sing_handler,
true);
6884 (*current_liboctave_error_handler) (
"unknown matrix type");
6889 #if defined (USE_QRSOLVE)
6890 retval =
qrsolve (*
this, b, err);
6906 return solve (mattype, b, info, rcond,
nullptr);
6914 return solve (mattype, b, info, rcond,
nullptr);
6921 return solve (mattype, b, info, rcond,
nullptr);
6927 solve_singularity_handler sing_handler,
6928 bool singular_fallback)
const
6931 int typ = mattype.
type (
false);
6934 typ = mattype.
type (*
this);
6937 retval = dsolve (mattype, b, err, rcond, sing_handler,
false);
6939 retval = utsolve (mattype, b, err, rcond, sing_handler,
false);
6941 retval = ltsolve (mattype, b, err, rcond, sing_handler,
false);
6943 retval = bsolve (mattype, b, err, rcond, sing_handler,
false);
6946 retval = trisolve (mattype, b, err, rcond, sing_handler,
false);
6948 retval = fsolve (mattype, b, err, rcond, sing_handler,
true);
6950 (*current_liboctave_error_handler) (
"unknown matrix type");
6955 #if defined (USE_QRSOLVE)
6956 retval =
qrsolve (*
this, b, err);
6970 return solve (mattype, b, info, rcond);
6978 return solve (mattype, b, info, rcond);
6985 return solve (mattype, b, info, rcond,
nullptr);
6991 solve_singularity_handler sing_handler)
const
6994 return solve (mattype, tmp, info, rcond,
7004 return solve (mattype, b, info, rcond,
nullptr);
7012 return solve (mattype, b, info, rcond,
nullptr);
7019 return solve (mattype, b, info, rcond,
nullptr);
7025 solve_singularity_handler sing_handler)
const
7028 return solve (mattype, tmp, info, rcond,
7037 return solve (b, info, rcond,
nullptr);
7044 return solve (b, info, rcond,
nullptr);
7049 double& rcond)
const
7051 return solve (b, info, rcond,
nullptr);
7057 solve_singularity_handler sing_handler)
const
7060 return solve (mattype, b, err, rcond, sing_handler);
7068 return solve (b, info, rcond,
nullptr);
7076 return solve (b, info, rcond,
nullptr);
7083 return solve (b, info, rcond,
nullptr);
7089 solve_singularity_handler sing_handler)
const
7092 return solve (mattype, b, err, rcond, sing_handler);
7100 return solve (b, info, rcond,
nullptr);
7107 return solve (b, info, rcond,
nullptr);
7113 solve_singularity_handler sing_handler)
const
7116 return solve (mattype, b, err, rcond, sing_handler);
7124 return solve (b, info, rcond,
nullptr);
7132 return solve (b, info, rcond,
nullptr);
7139 return solve (b, info, rcond,
nullptr);
7145 solve_singularity_handler sing_handler)
const
7148 return solve (mattype, b, err, rcond, sing_handler);
7155 return solve (b, info, rcond);
7162 return solve (b, info, rcond);
7167 double& rcond)
const
7169 return solve (b, info, rcond,
nullptr);
7175 solve_singularity_handler sing_handler)
const
7178 return solve (tmp, info, rcond,
7187 return solve (b, info, rcond,
nullptr);
7195 return solve (b, info, rcond,
nullptr);
7200 double& rcond)
const
7202 return solve (b, info, rcond,
nullptr);
7208 solve_singularity_handler sing_handler)
const
7211 return solve (tmp, info, rcond,
7236 if (jj <
cidx (i+1) &&
ridx (jj) == j)
7333 double r_val = val.real ();
7334 double i_val = val.imag ();
7336 if (r_val > max_val)
7339 if (i_val > max_val)
7342 if (r_val < min_val)
7345 if (i_val < min_val)
7391 if ((
rows () == 1 && dim == -1) || dim == 1)
7396 (
cidx (j+1) -
cidx (j) < nr ? 0.0 : 1.0), 1.0);
7410 Complex d = data (i); \
7411 tmp[ridx (i)] += d * conj (d)
7414 Complex d = data (i); \
7415 tmp[j] += d * conj (d)
7437 retval.
data (i) = std::abs (
data (i));
7462 os << a.
ridx (i) + 1 <<
' ' << j + 1 <<
' ';
7463 octave::write_value<Complex> (os, a.
data (i));
7476 return read_sparse_matrix<elt_type> (is, a, octave::read_value<Complex>);
7561 return do_mul_dm_sm<SparseComplexMatrix> (
d, a);
7566 return do_mul_sm_dm<SparseComplexMatrix> (a,
d);
7572 return do_mul_dm_sm<SparseComplexMatrix> (
d, a);
7577 return do_mul_sm_dm<SparseComplexMatrix> (a,
d);
7583 return do_mul_dm_sm<SparseComplexMatrix> (
d, a);
7588 return do_mul_sm_dm<SparseComplexMatrix> (a,
d);
7594 return do_add_dm_sm<SparseComplexMatrix> (
d, a);
7599 return do_add_dm_sm<SparseComplexMatrix> (
d, a);
7604 return do_add_dm_sm<SparseComplexMatrix> (
d, a);
7609 return do_add_sm_dm<SparseComplexMatrix> (a,
d);
7614 return do_add_sm_dm<SparseComplexMatrix> (a,
d);
7619 return do_add_sm_dm<SparseComplexMatrix> (a,
d);
7625 return do_sub_dm_sm<SparseComplexMatrix> (
d, a);
7630 return do_sub_dm_sm<SparseComplexMatrix> (
d, a);
7635 return do_sub_dm_sm<SparseComplexMatrix> (
d, a);
7640 return do_sub_sm_dm<SparseComplexMatrix> (a,
d);
7645 return do_sub_sm_dm<SparseComplexMatrix> (a,
d);
7650 return do_sub_sm_dm<SparseComplexMatrix> (a,
d);
7669 #define EMPTY_RETURN_CHECK(T) \
7670 if (nr == 0 || nc == 0) \
7713 if (a_nr == b_nr && a_nc == b_nc)
7723 bool ja_lt_max = ja < ja_max;
7727 bool jb_lt_max = jb < jb_max;
7729 while (ja_lt_max || jb_lt_max)
7732 if ((! jb_lt_max) || (ja_lt_max && (a.
ridx (ja) < b.
ridx (jb))))
7737 r.ridx (jx) = a.
ridx (ja);
7742 ja_lt_max= ja < ja_max;
7744 else if ((! ja_lt_max)
7745 || (jb_lt_max && (b.
ridx (jb) < a.
ridx (ja))))
7750 r.ridx (jx) = b.
ridx (jb);
7755 jb_lt_max= jb < jb_max;
7763 r.ridx (jx) = a.
ridx (ja);
7767 ja_lt_max= ja < ja_max;
7769 jb_lt_max= jb < jb_max;
7775 r.maybe_compress ();
7779 if (a_nr == 0 || a_nc == 0)
7780 r.resize (a_nr, a_nc);
7781 else if (b_nr == 0 || b_nc == 0)
7782 r.resize (b_nr, b_nc);
7830 if (a_nr == b_nr && a_nc == b_nc)
7840 bool ja_lt_max = ja < ja_max;
7844 bool jb_lt_max = jb < jb_max;
7846 while (ja_lt_max || jb_lt_max)
7849 if ((! jb_lt_max) || (ja_lt_max && (a.
ridx (ja) < b.
ridx (jb))))
7854 r.ridx (jx) = a.
ridx (ja);
7859 ja_lt_max= ja < ja_max;
7861 else if ((! ja_lt_max)
7862 || (jb_lt_max && (b.
ridx (jb) < a.
ridx (ja))))
7867 r.ridx (jx) = b.
ridx (jb);
7872 jb_lt_max= jb < jb_max;
7880 r.ridx (jx) = a.
ridx (ja);
7884 ja_lt_max= ja < ja_max;
7886 jb_lt_max= jb < jb_max;
7892 r.maybe_compress ();
7896 if (a_nr == 0 || a_nc == 0)
7897 r.resize (a_nr, a_nc);
7898 else if (b_nr == 0 || b_nc == 0)
7899 r.resize (b_nr, b_nc);
SparseComplexMatrix max(const Complex &c, const SparseComplexMatrix &m)
std::istream & operator>>(std::istream &is, SparseComplexMatrix &a)
std::ostream & operator<<(std::ostream &os, const SparseComplexMatrix &a)
SparseComplexMatrix operator-(const ComplexDiagMatrix &d, const SparseMatrix &a)
SparseComplexMatrix conj(const SparseComplexMatrix &a)
ComplexMatrix mul_trans(const ComplexMatrix &m, const SparseComplexMatrix &a)
SparseComplexMatrix operator*(const SparseComplexMatrix &m, const SparseMatrix &a)
#define EMPTY_RETURN_CHECK(T)
ComplexMatrix mul_herm(const ComplexMatrix &m, const SparseComplexMatrix &a)
SparseComplexMatrix min(const Complex &c, const SparseComplexMatrix &m)
SparseComplexMatrix operator+(const ComplexDiagMatrix &d, const SparseMatrix &a)
ComplexMatrix trans_mul(const SparseComplexMatrix &m, const ComplexMatrix &a)
ComplexMatrix herm_mul(const SparseComplexMatrix &m, const ComplexMatrix &a)
base_det< Complex > ComplexDET
#define SPARSE_SSM_BOOL_OPS(S, M)
#define SPARSE_ANY_OP(DIM)
#define SPARSE_BASE_REDUCTION_OP(RET_TYPE, EL_TYPE, ROW_EXPR, COL_EXPR, INIT_VAL, MT_RESULT)
#define SPARSE_FULL_MUL(RET_TYPE, EL_TYPE)
#define FULL_SPARSE_MUL_TRANS(RET_TYPE, EL_TYPE, CONJ_OP)
#define SPARSE_SSM_CMP_OPS(S, M)
#define SPARSE_SMSM_CMP_OPS(M1, M2)
#define SPARSE_REDUCTION_OP(RET_TYPE, EL_TYPE, OP, INIT_VAL, MT_RESULT)
#define SPARSE_SMS_CMP_OPS(M, S)
#define SPARSE_SPARSE_MUL(RET_TYPE, RET_EL_TYPE, EL_TYPE)
#define SPARSE_CUMSUM(RET_TYPE, ELT_TYPE, FCN)
#define SPARSE_CUMPROD(RET_TYPE, ELT_TYPE, FCN)
#define SPARSE_FULL_TRANS_MUL(RET_TYPE, EL_TYPE, CONJ_OP)
#define FULL_SPARSE_MUL(RET_TYPE, EL_TYPE)
#define SPARSE_SMSM_BOOL_OPS(M1, M2)
#define SPARSE_SMS_BOOL_OPS(M, S)
#define SPARSE_ALL_OP(DIM)
SM octinternal_do_mul_pm_sm(const PermMatrix &p, const SM &a)
SM octinternal_do_mul_sm_pm(const SM &a, const PermMatrix &p)
T & elem(octave_idx_type n)
Size of the specified dimension.
T * fortran_vec()
Size of the specified dimension.
octave_idx_type rows() const
void resize(const dim_vector &dv, const T &rfv)
Size of the specified dimension.
const T * data() const
Size of the specified dimension.
octave_idx_type cols() const
T & xelem(octave_idx_type n)
Size of the specified dimension.
void resize(octave_idx_type nr, octave_idx_type nc, const Complex &rfv=Complex(0))
ComplexColumnVector column(octave_idx_type i) const
octave_idx_type length() const
octave_idx_type cols() const
MSparse< T > permute(const Array< octave_idx_type > &vec, bool inv=false) const
MSparse< T > squeeze() const
MSparse< T > ipermute(const Array< octave_idx_type > &vec) const
MSparse< T > reshape(const dim_vector &new_dims) const
MSparse< T > & insert(const Sparse< T > &a, octave_idx_type r, octave_idx_type c)
MSparse< T > diag(octave_idx_type k=0) const
void mark_as_unsymmetric()
octave_idx_type * triangular_perm() const
MatrixType transpose() const
int type(bool quiet=true)
void mark_as_rectangular()
SparseComplexMatrix min(int dim=-1) const
ComplexRowVector row(octave_idx_type i) const
bool too_large_for_float() const
bool any_element_is_inf_or_nan() const
friend SparseComplexMatrix conj(const SparseComplexMatrix &a)
SparseComplexMatrix hermitian() const
SparseComplexMatrix sumsq(int dim=-1) const
SparseComplexMatrix diag(octave_idx_type k=0) const
ComplexColumnVector column(octave_idx_type i) const
ComplexMatrix matrix_value() const
bool all_elements_are_real() const
SparseComplexMatrix & insert(const SparseComplexMatrix &a, octave_idx_type r, octave_idx_type c)
SparseComplexMatrix max(int dim=-1) const
bool any_element_is_nan() const
SparseComplexMatrix inverse() const
bool all_integers(double &max_val, double &min_val) const
SparseComplexMatrix reshape(const dim_vector &new_dims) const
SparseComplexMatrix permute(const Array< octave_idx_type > &vec, bool inv=false) const
bool operator!=(const SparseComplexMatrix &a) const
SparseBoolMatrix any(int dim=-1) const
SparseComplexMatrix prod(int dim=-1) const
SparseComplexMatrix sum(int dim=-1) const
SparseComplexMatrix transpose() const
SparseComplexMatrix ipermute(const Array< octave_idx_type > &vec) const
bool operator==(const SparseComplexMatrix &a) const
SparseComplexMatrix concat(const SparseComplexMatrix &rb, const Array< octave_idx_type > &ra_idx)
SparseComplexMatrix squeeze() const
ComplexDET determinant() const
SparseBoolMatrix operator!() const
ComplexMatrix solve(MatrixType &mattype, const Matrix &b) const
SparseComplexMatrix cumsum(int dim=-1) const
SparseComplexMatrix cumprod(int dim=-1) const
SparseBoolMatrix all(int dim=-1) const
octave_idx_type cols() const
octave_idx_type * xcidx()
Array< T > array_value() const
octave_idx_type columns() const
bool test_any(F fcn) const
Sparse< T, Alloc > maybe_compress(bool remove_zeros=false)
T & elem(octave_idx_type n)
octave_idx_type nnz() const
Actual number of nonzero terms.
octave_idx_type rows() const
void change_capacity(octave_idx_type nz)
octave_idx_type * xridx()
Vector representing the dimensions (size) of an Array.
octave_idx_type ndims() const
Number of dimensions.
int first_non_singleton(int def=0) const
ColumnVector real(const ComplexColumnVector &a)
#define F77_DBLE_CMPLX_ARG(x)
#define F77_XFCN(f, F, args)
octave_f77_int_type F77_INT
void err_nonconformant(const char *op, octave_idx_type op1_len, octave_idx_type op2_len)
void err_nan_to_logical_conversion()
void warn_singular_matrix(double rcond)
OCTAVE_NORETURN liboctave_error_handler current_liboctave_error_handler
F77_RET_T const F77_INT F77_CMPLX const F77_INT F77_CMPLX * B
F77_RET_T const F77_INT const F77_INT const F77_INT F77_DBLE const F77_INT F77_DBLE const F77_INT F77_DBLE * Q
F77_RET_T const F77_INT F77_CMPLX * A
F77_RET_T const F77_DBLE const F77_DBLE F77_DBLE * d
bool too_large_for_float(double x)
bool mx_inline_all_real(std::size_t n, const std::complex< T > *x)
std::complex< double > Complex
#define OCTAVE_LOCAL_BUFFER(T, buf, size)
suitesparse_integer * to_suitesparse_intptr(octave_idx_type *i)
#define UMFPACK_ZNAME(name)
#define CHOLMOD_NAME(name)
octave_idx_type from_suitesparse_long(SuiteSparse_long x)
octave_idx_type from_size_t(std::size_t x)
const octave_base_value const Array< octave_idx_type > & ra_idx
template int8_t abs(int8_t)
template ComplexMatrix dmsolve< ComplexMatrix, SparseComplexMatrix, ComplexMatrix >(const SparseComplexMatrix &, const ComplexMatrix &, octave_idx_type &)
template ComplexMatrix dmsolve< ComplexMatrix, SparseComplexMatrix, Matrix >(const SparseComplexMatrix &, const Matrix &, octave_idx_type &)
RT dmsolve(const ST &a, const T &b, octave_idx_type &info)
template SparseComplexMatrix dmsolve< SparseComplexMatrix, SparseComplexMatrix, SparseMatrix >(const SparseComplexMatrix &, const SparseMatrix &, octave_idx_type &)
Matrix qrsolve(const SparseMatrix &a, const MArray< double > &b, octave_idx_type &info)
void SparseCholError(int status, char *file, int line, char *message)
int SparseCholPrint(const char *fmt,...)