Mercurial > hg > octave-lyh

comparison liboctave/CMatrix.cc @ 4330:9f86c2055b58

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 2003-02-18 22:31:55 by jwe]
author jwe
date Tue, 18 Feb 2003 22:31:55 +0000
parents d53c33d93440
children a6c22c2c9b09
comparison
equal deleted inserted replaced
4329:d53c33d93440 4330:9f86c2055b58
958 958
959 retval = *this; 959 retval = *this;
960 Complex *tmp_data = retval.fortran_vec (); 960 Complex *tmp_data = retval.fortran_vec ();
961 961
962 Array<Complex> z(1); 962 Array<Complex> z(1);
963 int lwork = 1; 963 int lwork = -1;
964 964
965 // Query the optimum work array size 965 // Query the optimum work array size.
966 966
967 F77_XFCN (zgetri, ZGETRI, (nc, tmp_data, nr, pipvt, 967 F77_XFCN (zgetri, ZGETRI, (nc, tmp_data, nr, pipvt,
968 z.fortran_vec (), lwork, info)); 968 z.fortran_vec (), lwork, info));
969 969
970 if (f77_exception_encountered) 970 if (f77_exception_encountered)
979 z.resize (lwork); 979 z.resize (lwork);
980 Complex *pz = z.fortran_vec (); 980 Complex *pz = z.fortran_vec ();
981 981
982 info = 0; 982 info = 0;
983 983
984 /* Calculate the norm of the matrix, for later use */ 984 // Calculate the norm of the matrix, for later use.
985 double anorm; 985 double anorm;
986 if (calc_cond) 986 if (calc_cond)
987 anorm = retval.abs().sum().row(0).max(); 987 anorm = retval.abs().sum().row(0).max();
988 988
989 F77_XFCN (zgetrf, ZGETRF, (nc, nc, tmp_data, nr, pipvt, info)); 989 F77_XFCN (zgetrf, ZGETRF, (nc, nc, tmp_data, nr, pipvt, info));
990 990
991 if (f77_exception_encountered) 991 if (f77_exception_encountered)
992 (*current_liboctave_error_handler) ("unrecoverable error in zgetrf"); 992 (*current_liboctave_error_handler) ("unrecoverable error in zgetrf");
993 else 993 else
994 { 994 {
995 /* Throw-away extra info LAPACK gives so as to not change output */ 995 // Throw-away extra info LAPACK gives so as to not change output.
996 rcond = 0.; 996 rcond = 0.;
997 if ( info != 0) 997 if ( info != 0)
998 info = -1; 998 info = -1;
999 else if (calc_cond) 999 else if (calc_cond)
1000 { 1000 {
1001 /* Now calculate the condition number for non-singular matrix */ 1001 // Now calculate the condition number for non-singular matrix.
1002 char job = '1'; 1002 char job = '1';
1003 Array<double> rz (2 * nc); 1003 Array<double> rz (2 * nc);
1004 double *prz = rz.fortran_vec (); 1004 double *prz = rz.fortran_vec ();
1005 F77_XFCN (zgecon, ZGECON, ( &job, nc, tmp_data, nr, anorm, 1005 F77_XFCN (zgecon, ZGECON, ( &job, nc, tmp_data, nr, anorm,
1006 rcond, pz, prz, info)); 1006 rcond, pz, prz, info));
1434 ComplexMatrix atmp = *this; 1434 ComplexMatrix atmp = *this;
1435 Complex *tmp_data = atmp.fortran_vec (); 1435 Complex *tmp_data = atmp.fortran_vec ();
1436 1436
1437 info = 0; 1437 info = 0;
1438 1438
1439 /* Calculate the norm of the matrix, for later use */ 1439 // Calculate the norm of the matrix, for later use.
1440 double anorm = 0; 1440 double anorm = 0;
1441 if (calc_cond) 1441 if (calc_cond)
1442 anorm = atmp.abs().sum().row(0).max(); 1442 anorm = atmp.abs().sum().row(0).max();
1443 1443
1444 F77_XFCN (zgetrf, ZGETRF, (nr, nc, tmp_data, nr, pipvt, info)); 1444 F77_XFCN (zgetrf, ZGETRF, (nr, nc, tmp_data, nr, pipvt, info));
1445 1445
1446 if (f77_exception_encountered) 1446 if (f77_exception_encountered)
1447 (*current_liboctave_error_handler) ("unrecoverable error in zgetrf"); 1447 (*current_liboctave_error_handler) ("unrecoverable error in zgetrf");
1448 else 1448 else
1449 { 1449 {
1450 /* Throw-away extra info LAPACK gives so as to not change output */ 1450 // Throw-away extra info LAPACK gives so as to not change output.
1451 rcond = 0.; 1451 rcond = 0.;
1452 if ( info != 0) 1452 if ( info != 0)
1453 { 1453 {
1454 info = -1; 1454 info = -1;
1455 retval = ComplexDET (); 1455 retval = ComplexDET ();
1456 } 1456 }
1457 else 1457 else
1458 { 1458 {
1459 if (calc_cond) 1459 if (calc_cond)
1460 { 1460 {
1461 /* Now calc the condition number for non-singular matrix */ 1461 // Now calc the condition number for non-singular matrix.
1462 char job = '1'; 1462 char job = '1';
1463 Array<Complex> z (2*nr); 1463 Array<Complex> z (2*nr);
1464 Complex *pz = z.fortran_vec (); 1464 Complex *pz = z.fortran_vec ();
1465 Array<double> rz (2*nr); 1465 Array<double> rz (2*nr);
1466 double *prz = rz.fortran_vec (); 1466 double *prz = rz.fortran_vec ();
1581 Array<Complex> z (2 * nc); 1581 Array<Complex> z (2 * nc);
1582 Complex *pz = z.fortran_vec (); 1582 Complex *pz = z.fortran_vec ();
1583 Array<double> rz (2 * nc); 1583 Array<double> rz (2 * nc);
1584 double *prz = rz.fortran_vec (); 1584 double *prz = rz.fortran_vec ();
1585 1585
1586 /* Calculate the norm of the matrix, for later use */ 1586 // Calculate the norm of the matrix, for later use.
1587 double anorm = atmp.abs().sum().row(0).max(); 1587 double anorm = atmp.abs().sum().row(0).max();
1588 1588
1589 F77_XFCN (zgetrf, ZGETRF, (nr, nr, tmp_data, nr, pipvt, info)); 1589 F77_XFCN (zgetrf, ZGETRF, (nr, nr, tmp_data, nr, pipvt, info));
1590 1590
1591 if (f77_exception_encountered) 1591 if (f77_exception_encountered)
1592 (*current_liboctave_error_handler) ("unrecoverable error in zgetrf"); 1592 (*current_liboctave_error_handler) ("unrecoverable error in zgetrf");
1593 else 1593 else
1594 { 1594 {
1595 /* Throw-away extra info LAPACK gives so as to not change output */ 1595 // Throw-away extra info LAPACK gives so as to not change output.
1596 rcond = 0.; 1596 rcond = 0.;
1597 if ( info != 0) 1597 if ( info != 0)
1598 { 1598 {
1599 info = -2; 1599 info = -2;
1600 1600
1605 ("matrix singular to machine precision"); 1605 ("matrix singular to machine precision");
1606 1606
1607 } 1607 }
1608 else 1608 else
1609 { 1609 {
1610 /* Now calculate the condition number for non-singular matrix */ 1610 // Now calculate the condition number for non-singular matrix.
1611 char job = '1'; 1611 char job = '1';
1612 F77_XFCN (zgecon, ZGECON, ( &job, nc, tmp_data, nr, anorm, 1612 F77_XFCN (zgecon, ZGECON, ( &job, nc, tmp_data, nr, anorm,
1613 rcond, pz, prz, info)); 1613 rcond, pz, prz, info));
1614 1614
1615 if (f77_exception_encountered) 1615 if (f77_exception_encountered)
1730 Array<Complex> z (2 * nc); 1730 Array<Complex> z (2 * nc);
1731 Complex *pz = z.fortran_vec (); 1731 Complex *pz = z.fortran_vec ();
1732 Array<double> rz (2 * nc); 1732 Array<double> rz (2 * nc);
1733 double *prz = rz.fortran_vec (); 1733 double *prz = rz.fortran_vec ();
1734 1734
1735 /* Calculate the norm of the matrix, for later use */ 1735 // Calculate the norm of the matrix, for later use.
1736 double anorm = atmp.abs().sum().row(0).max(); 1736 double anorm = atmp.abs().sum().row(0).max();
1737 1737
1738 F77_XFCN (zgetrf, ZGETRF, (nr, nr, tmp_data, nr, pipvt, info)); 1738 F77_XFCN (zgetrf, ZGETRF, (nr, nr, tmp_data, nr, pipvt, info));
1739 1739
1740 if (f77_exception_encountered) 1740 if (f77_exception_encountered)
1741 (*current_liboctave_error_handler) ("unrecoverable error in zgetrf"); 1741 (*current_liboctave_error_handler) ("unrecoverable error in zgetrf");
1742 else 1742 else
1743 { 1743 {
1744 /* Throw-away extra info LAPACK gives so as to not change output */ 1744 // Throw-away extra info LAPACK gives so as to not change output.
1745 rcond = 0.; 1745 rcond = 0.;
1746 if ( info != 0) 1746 if ( info != 0)
1747 { 1747 {
1748 info = -2; 1748 info = -2;
1749 1749
1754 ("matrix singular to machine precision, rcond = %g", 1754 ("matrix singular to machine precision, rcond = %g",
1755 rcond); 1755 rcond);
1756 } 1756 }
1757 else 1757 else
1758 { 1758 {
1759 /* Now calculate the condition number for non-singular matrix */ 1759 // Now calculate the condition number for non-singular matrix.
1760 char job = '1'; 1760 char job = '1';
1761 F77_XFCN (zgecon, ZGECON, ( &job, nc, tmp_data, nr, anorm, 1761 F77_XFCN (zgecon, ZGECON, ( &job, nc, tmp_data, nr, anorm,
1762 rcond, pz, prz, info)); 1762 rcond, pz, prz, info));
1763 1763
1764 if (f77_exception_encountered) 1764 if (f77_exception_encountered)