25 #if defined (HAVE_CONFIG_H) 116 elem (
i, j) =
static_cast<unsigned char> (
a.elem (
i, j));
122 if (
rows () !=
a.rows () ||
cols () !=
a.cols ())
131 return !(*
this ==
a);
162 if (r < 0 || r >=
rows () || c < 0 || c + a_len >
cols ())
163 (*current_liboctave_error_handler) (
"range error for insert");
181 if (r < 0 || r + a_len >
rows () || c < 0 || c >=
cols ())
182 (*current_liboctave_error_handler) (
"range error for insert");
201 if (r < 0 || r + a_nr >
rows () || c < 0 || c + a_nc >
cols ())
202 (*current_liboctave_error_handler) (
"range error for insert");
225 if (nr > 0 && nc > 0)
244 if (r1 < 0 || r2 < 0 || c1 < 0 || c2 < 0
245 || r1 >= nr || r2 >= nr || c1 >= nc || c2 >= nc)
246 (*current_liboctave_error_handler) (
"range error for fill");
251 if (r2 >= r1 && c2 >= c1)
273 retval.insert (*
this, 0, 0);
274 retval.insert (
a, 0, nc_insert);
284 (*current_liboctave_error_handler) (
"row dimension mismatch for append");
288 retval.insert (*
this, 0, 0);
289 retval.insert (
a, 0, nc_insert);
298 if (nr !=
a.numel ())
303 retval.insert (*
this, 0, 0);
304 retval.insert (
a, 0, nc_insert);
318 retval.insert (*
this, 0, 0);
319 retval.insert (
a, 0, nc_insert);
333 retval.insert (*
this, 0, 0);
334 retval.insert (
a, nr_insert, 0);
343 if (nc !=
a.numel ())
348 retval.insert (*
this, 0, 0);
349 retval.insert (
a, nr_insert, 0);
359 (*current_liboctave_error_handler) (
"column dimension mismatch for stack");
363 retval.insert (*
this, 0, 0);
364 retval.insert (
a, nr_insert, 0);
378 retval.insert (*
this, 0, 0);
379 retval.insert (
a, nr_insert, 0);
436 double sum = colsum.
xelem (
i);
455 return inverse (mattype, info, rcon, 0, 0);
463 return inverse (mattype, info, rcon, 0, 0);
468 bool calc_cond)
const 471 return inverse (mattype, info, rcon, force, calc_cond);
479 return inverse (mattype, info, rcon, 0, 0);
486 return inverse (mattype, info, rcon, 0, 0);
491 bool force,
bool calc_cond)
const 498 if (nr != nc || nr == 0 || nc == 0)
499 (*current_liboctave_error_handler) (
"inverse requires square matrix");
501 int typ = mattype.
type ();
505 double *tmp_data =
retval.fortran_vec ();
509 F77_XFCN (dtrtri, DTRTRI, (F77_CONST_CHAR_ARG2 (&uplo, 1),
510 F77_CONST_CHAR_ARG2 (&udiag, 1),
511 nr, tmp_data, nr, tmp_info
513 F77_CHAR_ARG_LEN (1)));
529 F77_XFCN (dtrcon, DTRCON, (F77_CONST_CHAR_ARG2 (&job, 1),
530 F77_CONST_CHAR_ARG2 (&uplo, 1),
531 F77_CONST_CHAR_ARG2 (&udiag, 1),
532 nr, tmp_data, nr, rcon,
533 work, iwork, dtrcon_info
536 F77_CHAR_ARG_LEN (1)));
538 if (dtrcon_info != 0)
542 if (info == -1 && ! force)
550 bool force,
bool calc_cond)
const 557 if (nr != nc || nr == 0 || nc == 0)
558 (*current_liboctave_error_handler) (
"inverse requires square matrix");
561 F77_INT *pipvt = ipvt.fortran_vec ();
564 double *tmp_data =
retval.fortran_vec ();
572 F77_XFCN (dgetri, DGETRI, (nc, tmp_data, nr, pipvt,
577 lwork =
static_cast<F77_INT> (z(0));
578 lwork = (lwork < 4 * nc ? 4 * nc : lwork);
590 F77_XFCN (dgetrf, DGETRF, (nc, nc, tmp_data, nr, pipvt, tmp_info));
606 F77_XFCN (dgecon, DGECON, (F77_CONST_CHAR_ARG2 (&job, 1),
607 nc, tmp_data, nr, anorm,
608 rcon, pz, piz, dgecon_info
609 F77_CHAR_ARG_LEN (1)));
611 if (dgecon_info != 0)
615 if (info == -1 && ! force)
621 F77_XFCN (dgetri, DGETRI, (nc, tmp_data, nr, pipvt,
622 pz, lwork, dgetri_info));
624 if (dgetri_info != 0)
636 bool force,
bool calc_cond)
const 638 int typ = mattype.
type (
false);
642 typ = mattype.
type (*
this);
645 ret =
tinverse (mattype, info, rcon, force, calc_cond);
654 rcon =
chol.rcond ();
657 ret =
chol.inverse ();
664 ret =
finverse (mattype, info, rcon, force, calc_cond);
666 if ((calc_cond || mattype.
ishermitian ()) && rcon == 0.
667 && ! (
numel () == 1))
694 * std::numeric_limits<double>::epsilon ();
700 while (r >= 0 && sigma.
elem (r) < tol)
704 return Matrix (nc, nr, 0.0);
709 Matrix Vr =
V.extract (0, 0, nc-1, r);
714 #if defined (HAVE_FFTW) 724 size_t npts, nsamples;
726 if (nr == 1 || nc == 1)
728 npts = (nr > nc ? nr : nc);
753 size_t npts, nsamples;
755 if (nr == 1 || nc == 1)
757 npts = (nr > nc ? nr : nc);
815 if (nr == 1 || nc == 1)
817 npts = octave::to_f77_int (nr > nc ? nr : nc);
822 npts = octave::to_f77_int (nr);
829 Complex *pwsave = wsave.fortran_vec ();
858 if (nr == 1 || nc == 1)
860 npts = octave::to_f77_int (nr > nc ? nr : nc);
865 npts = octave::to_f77_int (nr);
872 Complex *pwsave = wsave.fortran_vec ();
888 tmp_data[j] = tmp_data[j] / static_cast<double> (npts);
903 if (nr == 1 || nc == 1)
905 npts = (nr > nc ? nr : nc);
917 Complex *pwsave = wsave.fortran_vec ();
924 for (
F77_INT j = 0; j < nsamples; j++)
937 pwsave = wsave.fortran_vec ();
944 for (
F77_INT j = 0; j < nsamples; j++)
949 prow[
i] = tmp_data[
i*nr + j];
955 tmp_data[
i*nr + j] = prow[
i];
971 if (nr == 1 || nc == 1)
973 npts = (nr > nc ? nr : nc);
985 Complex *pwsave = wsave.fortran_vec ();
992 for (
F77_INT j = 0; j < nsamples; j++)
1000 for (
F77_INT j = 0; j < npts*nsamples; j++)
1001 tmp_data[j] = tmp_data[j] / static_cast<double> (npts);
1008 pwsave = wsave.fortran_vec ();
1015 for (
F77_INT j = 0; j < nsamples; j++)
1020 prow[
i] = tmp_data[
i*nr + j];
1026 tmp_data[
i*nr + j] = prow[
i] / static_cast<double> (npts);
1053 return determinant (mattype, info, rcon, calc_cond);
1069 (*current_liboctave_error_handler) (
"matrix must be square");
1071 volatile int typ = mattype.
type ();
1077 typ = mattype.
type (*
this);
1093 anorm =
norm1 (*
this);
1098 F77_XFCN (dpotrf, DPOTRF, (F77_CONST_CHAR_ARG2 (&job, 1), nr,
1099 tmp_data, nr, tmp_info
1100 F77_CHAR_ARG_LEN (1)));
1117 F77_XFCN (dpocon, DPOCON, (F77_CONST_CHAR_ARG2 (&job, 1),
1118 nr, tmp_data, nr, anorm,
1119 rcon, pz, piz, tmp_info
1120 F77_CHAR_ARG_LEN (1)));
1134 (*current_liboctave_error_handler) (
"det: invalid dense matrix type");
1150 anorm =
norm1 (*
this);
1152 F77_XFCN (dgetrf, DGETRF, (nr, nr, tmp_data, nr, pipvt, tmp_info));
1174 F77_XFCN (dgecon, DGECON, (F77_CONST_CHAR_ARG2 (&job, 1),
1175 nc, tmp_data, nr, anorm,
1176 rcon, pz, piz, tmp_info
1177 F77_CHAR_ARG_LEN (1)));
1191 double c = atmp(
i,
i);
1205 return rcond (mattype);
1216 (*current_liboctave_error_handler) (
"matrix must be square");
1218 if (nr == 0 || nc == 0)
1222 volatile int typ = mattype.
type ();
1225 typ = mattype.
type (*
this);
1241 F77_XFCN (dtrcon, DTRCON, (F77_CONST_CHAR_ARG2 (&
norm, 1),
1242 F77_CONST_CHAR_ARG2 (&uplo, 1),
1243 F77_CONST_CHAR_ARG2 (&dia, 1),
1244 nr, tmp_data, nr, rcon,
1246 F77_CHAR_ARG_LEN (1)
1247 F77_CHAR_ARG_LEN (1)
1248 F77_CHAR_ARG_LEN (1)));
1254 (*current_liboctave_error_handler)
1255 (
"permuted triangular matrix not implemented");
1269 F77_XFCN (dtrcon, DTRCON, (F77_CONST_CHAR_ARG2 (&
norm, 1),
1270 F77_CONST_CHAR_ARG2 (&uplo, 1),
1271 F77_CONST_CHAR_ARG2 (&dia, 1),
1272 nr, tmp_data, nr, rcon,
1274 F77_CHAR_ARG_LEN (1)
1275 F77_CHAR_ARG_LEN (1)
1276 F77_CHAR_ARG_LEN (1)));
1282 (*current_liboctave_error_handler)
1283 (
"permuted triangular matrix not implemented");
1286 double anorm = -1.0;
1296 anorm =
norm1 (atmp);
1298 F77_XFCN (dpotrf, DPOTRF, (F77_CONST_CHAR_ARG2 (&job, 1), nr,
1300 F77_CHAR_ARG_LEN (1)));
1315 F77_XFCN (dpocon, DPOCON, (F77_CONST_CHAR_ARG2 (&job, 1),
1316 nr, tmp_data, nr, anorm,
1318 F77_CHAR_ARG_LEN (1)));
1336 anorm =
norm1 (atmp);
1343 F77_XFCN (dgetrf, DGETRF, (nr, nr, tmp_data, nr, pipvt, info));
1353 F77_XFCN (dgecon, DGECON, (F77_CONST_CHAR_ARG2 (&job, 1),
1354 nc, tmp_data, nr, anorm,
1356 F77_CHAR_ARG_LEN (1)));
1372 double& rcon, solve_singularity_handler sing_handler,
1384 (*current_liboctave_error_handler)
1385 (
"matrix dimension mismatch solution of linear equations");
1387 if (nr == 0 || nc == 0 ||
b_nc == 0)
1391 volatile int typ = mattype.
type ();
1394 (*current_liboctave_error_handler) (
"incorrect matrix type");
1400 (*current_liboctave_error_handler)
1401 (
"permuted triangular matrix not implemented");
1414 F77_XFCN (dtrtrs, DTRTRS, (F77_CONST_CHAR_ARG2 (&uplo, 1),
1415 F77_CONST_CHAR_ARG2 (&trans, 1),
1416 F77_CONST_CHAR_ARG2 (&dia, 1),
1417 nr,
b_nc, tmp_data, nr,
1419 F77_CHAR_ARG_LEN (1)
1420 F77_CHAR_ARG_LEN (1)
1421 F77_CHAR_ARG_LEN (1)));
1436 F77_XFCN (dtrcon, DTRCON, (F77_CONST_CHAR_ARG2 (&
norm, 1),
1437 F77_CONST_CHAR_ARG2 (&uplo, 1),
1438 F77_CONST_CHAR_ARG2 (&dia, 1),
1439 nr, tmp_data, nr, rcon,
1441 F77_CHAR_ARG_LEN (1)
1442 F77_CHAR_ARG_LEN (1)
1443 F77_CHAR_ARG_LEN (1)));
1450 volatile double rcond_plus_one = rcon + 1.0;
1457 sing_handler (rcon);
1469 double& rcon, solve_singularity_handler sing_handler,
1481 (*current_liboctave_error_handler)
1482 (
"matrix dimension mismatch solution of linear equations");
1484 if (nr == 0 || nc == 0 ||
b_nc == 0)
1488 volatile int typ = mattype.
type ();
1491 (*current_liboctave_error_handler) (
"incorrect matrix type");
1497 (*current_liboctave_error_handler)
1498 (
"permuted triangular matrix not implemented");
1511 F77_XFCN (dtrtrs, DTRTRS, (F77_CONST_CHAR_ARG2 (&uplo, 1),
1512 F77_CONST_CHAR_ARG2 (&trans, 1),
1513 F77_CONST_CHAR_ARG2 (&dia, 1),
1514 nr,
b_nc, tmp_data, nr,
1516 F77_CHAR_ARG_LEN (1)
1517 F77_CHAR_ARG_LEN (1)
1518 F77_CHAR_ARG_LEN (1)));
1533 F77_XFCN (dtrcon, DTRCON, (F77_CONST_CHAR_ARG2 (&
norm, 1),
1534 F77_CONST_CHAR_ARG2 (&uplo, 1),
1535 F77_CONST_CHAR_ARG2 (&dia, 1),
1536 nr, tmp_data, nr, rcon,
1538 F77_CHAR_ARG_LEN (1)
1539 F77_CHAR_ARG_LEN (1)
1540 F77_CHAR_ARG_LEN (1)));
1547 volatile double rcond_plus_one = rcon + 1.0;
1554 sing_handler (rcon);
1566 double& rcon, solve_singularity_handler sing_handler,
1567 bool calc_cond)
const 1574 if (nr != nc || nr !=
b.rows ())
1576 (
"matrix dimension mismatch solution of linear equations");
1578 if (nr == 0 ||
b.cols () == 0)
1582 volatile int typ = mattype.
type ();
1595 anorm =
norm1 (atmp);
1599 F77_XFCN (dpotrf, DPOTRF, (F77_CONST_CHAR_ARG2 (&job, 1), nr,
1600 tmp_data, nr, tmp_info
1601 F77_CHAR_ARG_LEN (1)));
1623 F77_XFCN (dpocon, DPOCON, (F77_CONST_CHAR_ARG2 (&job, 1),
1624 nr, tmp_data, nr, anorm,
1625 rcon, pz, piz, tmp_info
1626 F77_CHAR_ARG_LEN (1)));
1633 volatile double rcond_plus_one = rcon + 1.0;
1640 sing_handler (rcon);
1654 F77_XFCN (dpotrs, DPOTRS, (F77_CONST_CHAR_ARG2 (&job, 1),
1655 nr,
b_nc, tmp_data, nr,
1657 F77_CHAR_ARG_LEN (1)));
1680 anorm =
norm1 (atmp);
1689 F77_XFCN (dgetrf, DGETRF, (nr, nr, tmp_data, nr, pipvt, tmp_info));
1700 sing_handler (rcon);
1713 F77_XFCN (dgecon, DGECON, (F77_CONST_CHAR_ARG2 (&job, 1),
1714 nc, tmp_data, nr, anorm,
1715 rcon, pz, piz, tmp_info
1716 F77_CHAR_ARG_LEN (1)));
1723 volatile double rcond_plus_one = rcon + 1.0;
1730 sing_handler (rcon);
1745 F77_XFCN (dgetrs, DGETRS, (F77_CONST_CHAR_ARG2 (&job, 1),
1746 nr,
b_nc, tmp_data, nr,
1748 F77_CHAR_ARG_LEN (1)));
1757 (*current_liboctave_error_handler) (
"incorrect matrix type");
1768 return solve (mattype,
b, info, rcon,
nullptr);
1776 return solve (mattype,
b, info, rcon,
nullptr);
1783 return solve (mattype,
b, info, rcon,
nullptr);
1788 double& rcon, solve_singularity_handler sing_handler,
1792 int typ = mattype.
type ();
1795 typ = mattype.
type (*
this);
1799 retval =
utsolve (mattype,
b, info, rcon, sing_handler,
true, transt);
1801 retval =
ltsolve (mattype,
b, info, rcon, sing_handler,
true, transt);
1806 retval =
fsolve (mattype,
b, info, rcon, sing_handler,
true);
1808 (*current_liboctave_error_handler) (
"unknown matrix type");
1825 return solve (mattype,
b, info, rcon,
nullptr);
1833 return solve (mattype,
b, info, rcon,
nullptr);
1840 return solve (mattype,
b, info, rcon,
nullptr);
1851 double *rd =
retval.fortran_vec ();
1867 const double *smd = sm.
data ();
1877 solve_singularity_handler sing_handler,
1881 tmp =
solve (mattype,
tmp, info, rcon, sing_handler, singular_fallback,
1890 return solve (mattype,
b, info, rcon);
1898 return solve (mattype,
b, info, rcon);
1905 return solve (mattype,
b, info, rcon,
nullptr);
1911 solve_singularity_handler sing_handler,
1915 tmp =
solve (mattype,
tmp, info, rcon, sing_handler,
true, transt);
1916 return tmp.column (static_cast<octave_idx_type> (0));
1923 return tmp.solve (mattype,
b);
1931 return tmp.solve (mattype,
b, info);
1939 return tmp.solve (mattype,
b, info, rcon);
1945 solve_singularity_handler sing_handler,
1949 return tmp.solve (mattype,
b, info, rcon, sing_handler, transt);
1957 return solve (
b, info, rcon,
nullptr);
1964 return solve (
b, info, rcon,
nullptr);
1970 return solve (
b, info, rcon,
nullptr);
1975 double& rcon, solve_singularity_handler sing_handler,
1979 return solve (mattype,
b, info, rcon, sing_handler,
true, transt);
1986 return tmp.solve (
b);
1993 return tmp.solve (
b, info);
2001 return tmp.solve (
b, info, rcon);
2006 solve_singularity_handler sing_handler,
2010 return tmp.solve (
b, info, rcon, sing_handler, transt);
2017 return solve (
b, info, rcon);
2024 return solve (
b, info, rcon);
2030 return solve (
b, info, rcon,
nullptr);
2035 solve_singularity_handler sing_handler,
2039 return solve (mattype,
b, info, rcon, sing_handler, transt);
2046 return tmp.solve (
b);
2053 return tmp.solve (
b, info);
2061 return tmp.solve (
b, info, rcon);
2067 solve_singularity_handler sing_handler,
2071 return tmp.solve (
b, info, rcon, sing_handler, transt);
2080 return lssolve (
b, info, rank, rcon);
2088 return lssolve (
b, info, rank, rcon);
2096 return lssolve (
b, info, rank, rcon);
2105 F77_INT nrhs = octave::to_f77_int (
b.cols ());
2114 (*current_liboctave_error_handler)
2115 (
"matrix dimension mismatch solution of linear equations");
2117 if (m == 0 || n == 0 ||
b_nc == 0)
2121 volatile F77_INT minmn = (m < n ? m : n);
2122 F77_INT maxmn = (m > n ? m : n);
2128 for (
F77_INT j = 0; j < nrhs; j++)
2138 double *pretval =
retval.fortran_vec ();
2140 double *ps =
s.fortran_vec ();
2149 F77_CONST_CHAR_ARG2 (
" ", 1),
2151 F77_CHAR_ARG_LEN (6)
2152 F77_CHAR_ARG_LEN (1));
2156 F77_CONST_CHAR_ARG2 (
" ", 1),
2157 m, n, nrhs, -1, mnthr
2158 F77_CHAR_ARG_LEN (6)
2159 F77_CHAR_ARG_LEN (1));
2163 double dminmn =
static_cast<double> (minmn);
2164 double dsmlsizp1 =
static_cast<double> (smlsiz+1);
2172 F77_INT liwork = 3 * minmn * nlvl + 11 * minmn;
2181 F77_XFCN (dgelsd, DGELSD, (m, n, nrhs, tmp_data, m, pretval, maxmn,
2183 lwork, piwork, tmp_info));
2192 if (n > m && n >= mnthr)
2195 = 9*m + 2*m*smlsiz + 8*m*nlvl + m*nrhs + (smlsiz+1)*(smlsiz+1);
2208 if (wlalsd > addend)
2211 const F77_INT lworkaround = 4*m + m*m + addend;
2213 if (work(0) < lworkaround)
2214 work(0) = lworkaround;
2219 = 12*n + 2*n*smlsiz + 8*n*nlvl + n*nrhs + (smlsiz+1)*(smlsiz+1);
2221 if (work(0) < lworkaround)
2222 work(0) = lworkaround;
2225 lwork =
static_cast<F77_INT> (work(0));
2228 anorm =
norm1 (*
this);
2243 F77_XFCN (dgelsd, DGELSD, (m, n, nrhs, tmp_data, m, pretval,
2244 maxmn, ps, rcon, tmp_rank,
2251 if (
s.elem (0) == 0.0)
2254 rcon =
s.elem (minmn - 1) /
s.elem (0);
2270 return tmp.lssolve (
b, info, rank, rcon);
2279 return tmp.lssolve (
b, info, rank, rcon);
2288 return tmp.lssolve (
b, info, rank, rcon);
2296 return tmp.lssolve (
b, info, rank, rcon);
2305 return lssolve (
b, info, rank, rcon);
2313 return lssolve (
b, info, rank, rcon);
2321 return lssolve (
b, info, rank, rcon);
2335 F77_INT b_nel = octave::to_f77_int (
b.numel ());
2338 (*current_liboctave_error_handler)
2339 (
"matrix dimension mismatch solution of linear equations");
2341 if (m == 0 || n == 0)
2345 volatile F77_INT minmn = (m < n ? m : n);
2346 F77_INT maxmn = (m > n ? m : n);
2362 double *pretval =
retval.fortran_vec ();
2364 double *ps =
s.fortran_vec ();
2373 F77_CONST_CHAR_ARG2 (
" ", 1),
2375 F77_CHAR_ARG_LEN (6)
2376 F77_CHAR_ARG_LEN (1));
2380 double dminmn =
static_cast<double> (minmn);
2381 double dsmlsizp1 =
static_cast<double> (smlsiz+1);
2388 F77_INT liwork = 3 * minmn * nlvl + 11 * minmn;
2397 F77_XFCN (dgelsd, DGELSD, (m, n, nrhs, tmp_data, m, pretval, maxmn,
2399 lwork, piwork, tmp_info));
2404 lwork =
static_cast<F77_INT> (work(0));
2407 F77_XFCN (dgelsd, DGELSD, (m, n, nrhs, tmp_data, m, pretval,
2408 maxmn, ps, rcon, tmp_rank,
2417 if (
s.elem (0) == 0.0)
2420 rcon =
s.elem (minmn - 1) /
s.elem (0);
2436 return tmp.lssolve (
b, info, rank, rcon);
2445 return tmp.lssolve (
b, info, rank, rcon);
2454 return tmp.lssolve (
b, info, rank, rcon);
2462 return tmp.lssolve (
b, info, rank, rcon);
2514 F77_INT a_len = octave::to_f77_int (
a.numel ());
2517 double *
c =
retval.fortran_vec ();
2519 F77_XFCN (dgemm, DGEMM, (F77_CONST_CHAR_ARG2 (
"N", 1),
2520 F77_CONST_CHAR_ARG2 (
"N", 1),
2521 len, a_len, 1, 1.0, v.
data (), len,
2522 a.data (), 1, 0.0,
c, len
2523 F77_CHAR_ARG_LEN (1)
2524 F77_CHAR_ARG_LEN (1)));
2596 if (nr == 1 || nc == 1)
2619 if (nr > 0 && nc > 0)
2630 for (idx_j = 0; idx_j < nc; idx_j++)
2632 tmp_min =
elem (
i, idx_j);
2644 else if (
tmp < tmp_min)
2674 if (nr > 0 && nc > 0)
2685 for (idx_j = 0; idx_j < nc; idx_j++)
2687 tmp_max =
elem (
i, idx_j);
2699 else if (
tmp > tmp_max)
2729 if (nr > 0 && nc > 0)
2740 for (idx_i = 0; idx_i < nr; idx_i++)
2742 tmp_min =
elem (idx_i, j);
2754 else if (
tmp < tmp_min)
2761 result.elem (j) = tmp_min;
2784 if (nr > 0 && nc > 0)
2795 for (idx_i = 0; idx_i < nr; idx_i++)
2797 tmp_max =
elem (idx_i, j);
2809 else if (
tmp > tmp_max)
2816 result.elem (j) = tmp_max;
2845 if (nr > 0 && nc > 0)
2851 tmp = octave_read_value<double> (
is);
2853 a.elem (
i, j) =
tmp;
2865 double cc,
s, temp_r;
2913 F77_XFCN (dtrsyl, DTRSYL, (F77_CONST_CHAR_ARG2 (
"N", 1),
2914 F77_CONST_CHAR_ARG2 (
"N", 1),
2917 F77_CHAR_ARG_LEN (1)
2918 F77_CHAR_ARG_LEN (1)));
2956 return trans ?
'T' :
'N';
2970 F77_INT a_nr = octave::to_f77_int (tra ?
a.cols () :
a.rows ());
2971 F77_INT a_nc = octave::to_f77_int (tra ?
a.rows () :
a.cols ());
2973 F77_INT b_nr = octave::to_f77_int (trb ?
b.cols () :
b.rows ());
2974 F77_INT b_nc = octave::to_f77_int (trb ?
b.rows () :
b.cols ());
2981 else if (
a.data () ==
b.data () &&
a_nr ==
b_nc && tra != trb)
2983 F77_INT lda = octave::to_f77_int (
a.rows ());
2986 double *
c =
retval.fortran_vec ();
2989 F77_XFCN (dsyrk, DSYRK, (F77_CONST_CHAR_ARG2 (
"U", 1),
2990 F77_CONST_CHAR_ARG2 (&ctra, 1),
2992 a.data (), lda, 0.0,
c,
a_nr 2993 F77_CHAR_ARG_LEN (1)
2994 F77_CHAR_ARG_LEN (1)));
2995 for (
int j = 0; j <
a_nr; j++)
2996 for (
int i = 0;
i < j;
i++)
3002 F77_INT lda = octave::to_f77_int (
a.rows ());
3003 F77_INT tda = octave::to_f77_int (
a.cols ());
3004 F77_INT ldb = octave::to_f77_int (
b.rows ());
3005 F77_INT tdb = octave::to_f77_int (
b.cols ());
3008 double *
c =
retval.fortran_vec ();
3017 F77_XFCN (dgemv, DGEMV, (F77_CONST_CHAR_ARG2 (&ctra, 1),
3018 lda, tda, 1.0,
a.data (), lda,
3019 b.data (), 1, 0.0,
c, 1
3020 F77_CHAR_ARG_LEN (1)));
3026 F77_XFCN (dgemv, DGEMV, (F77_CONST_CHAR_ARG2 (&crevtrb, 1),
3027 ldb, tdb, 1.0,
b.data (), ldb,
3028 a.data (), 1, 0.0,
c, 1
3029 F77_CHAR_ARG_LEN (1)));
3035 F77_XFCN (dgemm, DGEMM, (F77_CONST_CHAR_ARG2 (&ctra, 1),
3036 F77_CONST_CHAR_ARG2 (&ctrb, 1),
3038 lda,
b.data (), ldb, 0.0,
c,
a_nr 3039 F77_CHAR_ARG_LEN (1)
3040 F77_CHAR_ARG_LEN (1)));
3055 #define EMPTY_RETURN_CHECK(T) \ 3056 if (nr == 0 || nc == 0) \ 3105 if (nr !=
b.rows () || nc !=
b.columns ())
3107 (
"two-arg min requires same size arguments");
3169 if (nr !=
b.rows () || nc !=
b.columns ())
3171 (
"two-arg max requires same size arguments");
3194 if (x2.
numel () != m)
3196 (
"linspace: vectors must be of equal length");
3211 double *delta = &
retval(0, n-1);
3213 delta[
i] = (x2(
i) - x1(
i)) / (n - 1);
void octave_write_double(std::ostream &os, double d)
octave_idx_type rows(void) const
Matrix linspace(const ColumnVector &x1, const ColumnVector &x2, octave_idx_type n)
#define EMPTY_RETURN_CHECK(T)
Matrix pseudo_inverse(double tol=0.0) const
bool operator==(const Matrix &a) const
F77_RET_T const F77_INT const F77_INT const F77_INT const F77_DBLE const F77_DBLE F77_INT F77_DBLE * V
ComplexMatrix fourier2d(void) const
RowVector column_max(void) const
ColumnVector extract_diag(octave_idx_type k=0) const
Matrix inverse(void) const
#define MM_BOOL_OPS(M1, M2)
subroutine xddot(n, dx, incx, dy, incy, retval)
Matrix Givens(double x, double y)
static Matrix stack_complex_matrix(const ComplexMatrix &cm)
RowVector row(octave_idx_type i) const
const double * data(void) const
Matrix ltsolve(MatrixType &mattype, const Matrix &b, octave_idx_type &info, double &rcon, solve_singularity_handler sing_handler, bool calc_cond=false, blas_trans_type transt=blas_no_trans) const
static const idx_vector colon
DET determinant(void) const
subroutine zffti(n, wsave)
identity matrix If supplied two scalar respectively For allows like xample val
bool ishermitian(void) const
OCTAVE_NORETURN liboctave_error_handler current_liboctave_error_handler
#define F77_DBLE_CMPLX_ARG(x)
void mx_inline_real(size_t n, T *r, const std::complex< T > *x)
#define SM_BOOL_OPS(S, M)
ComplexMatrix fourier(void) const
const T * fortran_vec(void) const
friend class ComplexMatrix
Matrix Sylvester(const Matrix &a, const Matrix &b, const Matrix &c)
static char get_blas_trans_arg(bool trans)
Matrix sumsq(int dim=-1) const
Matrix min(double d, const Matrix &m)
static int ifft(const Complex *in, Complex *out, size_t npts, size_t nsamples=1, octave_idx_type stride=1, octave_idx_type dist=-1)
static int fftNd(const double *, Complex *, const int, const dim_vector &)
NDArray cumsum(int dim=-1) const
Matrix sum(int dim=-1) const
subroutine xilaenv(ispec, name, opts, n1, n2, n3, n4, retval)
octave_idx_type columns(void) const
Matrix stack(const Matrix &a) const
double & elem(octave_idx_type n)
bool operator!=(const Matrix &a) const
Matrix & operator+=(const DiagMatrix &a)
nd example oindent opens the file binary numeric values will be read assuming they are stored in IEEE format with the least significant bit and then converted to the native representation Opening a file that is already open simply opens it again and returns a separate file id It is not an error to open a file several though writing to the same file through several different file ids may produce unexpected results The possible values of text mode reading and writing automatically converts linefeeds to the appropriate line end character for the you may append a you must also open the file in binary mode The parameter conversions are currently only supported for and permissions will be set to and then everything is written in a single operation This is very efficient and improves performance c
#define F77_XFCN(f, F, args)
NDArray sumsq(int dim=-1) const
Complex log2(const Complex &x)
std::ostream & operator<<(std::ostream &os, const Matrix &a)
Matrix cumsum(int dim=-1) const
F77_RET_T const F77_REAL const F77_REAL F77_REAL &F77_RET_T const F77_DBLE const F77_DBLE F77_DBLE &F77_RET_T const F77_DBLE F77_DBLE &F77_RET_T const F77_REAL F77_REAL &F77_RET_T const F77_DBLE const F77_DBLE F77_DBLE * d
subroutine zfftb(n, c, wsave)
void mx_inline_imag(size_t n, T *r, const std::complex< T > *x)
octave_idx_type cols(void) const
Array< T > & insert(const Array< T > &a, const Array< octave_idx_type > &idx)
Insert an array into another at a specified position.
octave_value resize(const dim_vector &dv, bool fill=false) const
boolMatrix all(int dim=-1) const
calling an anonymous function involves an overhead quite comparable to the overhead of an m file function Passing a handle to a built in function is because the interpreter is not involved in the internal loop For a
boolNDArray all(int dim=-1) const
Matrix & fill(double val)
Matrix cumprod(int dim=-1) const
Matrix transpose(void) const
static int fft(const double *in, Complex *out, size_t npts, size_t nsamples=1, octave_idx_type stride=1, octave_idx_type dist=-1)
int type(bool quiet=true)
NDArray diag(octave_idx_type k=0) const
Matrix diag(octave_idx_type k=0) const
Array< double > index(const idx_vector &i) const
Indexing without resizing.
bool issquare(void) const
T schur_matrix(void) const
Matrix max(double d, const Matrix &m)
std::istream & operator>>(std::istream &is, Matrix &a)
Matrix utsolve(MatrixType &mattype, const Matrix &b, octave_idx_type &info, double &rcon, solve_singularity_handler sing_handler, bool calc_cond=false, blas_trans_type transt=blas_no_trans) const
double norm(const ColumnVector &v)
void resize(const dim_vector &dv, const T &rfv)
Resizing (with fill).
Matrix append(const Matrix &a) const
void err_nonconformant(const char *op, octave_idx_type op1_len, octave_idx_type op2_len)
ColumnVector row_max(void) const
RowVector column_min(void) const
void mark_as_unsymmetric(void)
T unitary_matrix(void) const
ColumnVector column(octave_idx_type i) const
Matrix xgemm(const Matrix &a, const Matrix &b, blas_trans_type transa, blas_trans_type transb)
ColumnVector extract(octave_idx_type r1, octave_idx_type r2) const
boolMatrix any(int dim=-1) const
Matrix real(const ComplexMatrix &a)
void mark_as_rectangular(void)
the exceeded dimensions are set to if fewer subscripts than dimensions are the exceeding dimensions are merged into the final requested dimension For consider the following dims
F77_RET_T F77_FUNC(xerbla, XERBLA)
Matrix & insert(const Matrix &a, octave_idx_type r, octave_idx_type c)
With real return the complex result
ComplexMatrix ifourier2d(void) const
NDArray cumprod(int dim=-1) const
static double norm1(const Matrix &a)
double & xelem(octave_idx_type n)
Matrix extract(octave_idx_type r1, octave_idx_type c1, octave_idx_type r2, octave_idx_type c2) const
Matrix & operator-=(const DiagMatrix &a)
static ComplexMatrix unstack_complex_matrix(const Matrix &sm)
Matrix tinverse(MatrixType &mattype, octave_idx_type &info, double &rcon, bool force, bool calc_cond) const
bool issymmetric(void) const
This is a simple wrapper template that will subclass an Array<T> type or any later type derived from ...
Matrix imag(const ComplexMatrix &a)
octave_f77_int_type F77_INT
OCTAVE_EXPORT octave_value_list or N dimensional array whose elements are all equal to the IEEE symbol NaN(Not a Number). NaN is the result of operations which do not produce a well defined 0 result. Common operations which produce a NaN are arithmetic with infinity ex($\infty - \infty$)
Matrix prod(int dim=-1) const
static int ifftNd(const Complex *, Complex *, const int, const dim_vector &)
Matrix solve(MatrixType &mattype, const Matrix &b) const
NDArray sum(int dim=-1) const
the element is set to zero In other the statement xample y
void scale(Matrix &m, double x, double y, double z)
ComplexMatrix ifourier(void) const
#define OCTAVE_LOCAL_BUFFER(T, buf, size)
#define MS_BOOL_OPS(M, S)
#define MM_CMP_OPS(M1, M2)
char get_blas_char(blas_trans_type transt)
NDArray prod(int dim=-1) const
boolNDArray any(int dim=-1) const
std::complex< double > Complex
Matrix extract_n(octave_idx_type r1, octave_idx_type c1, octave_idx_type nr, octave_idx_type nc) const
Matrix fsolve(MatrixType &mattype, const Matrix &b, octave_idx_type &info, double &rcon, solve_singularity_handler sing_handler, bool calc_cond=false) const
Matrix finverse(MatrixType &mattype, octave_idx_type &info, double &rcon, bool force, bool calc_cond) const
octave_idx_type numel(void) const
Number of elements in the array.
bool mx_inline_equal(size_t n, const T1 *x, const T2 *y)
write the output to stdout if nargout is
Vector representing the dimensions (size) of an Array.
void warn_singular_matrix(double rcond)
Matrix lssolve(const Matrix &b) const
F77_RET_T const F77_REAL const F77_REAL F77_REAL &F77_RET_T const F77_DBLE const F77_DBLE F77_DBLE &F77_RET_T const F77_DBLE F77_DBLE &F77_RET_T const F77_REAL F77_REAL &F77_RET_T const F77_DBLE * x
subroutine zfftf(n, c, wsave)
ColumnVector row_min(void) const
Matrix operator*(const ColumnVector &v, const RowVector &a)