# HG changeset patch # User jwe # Date 824094459 0 # Node ID 33392389484893ac04f62c26d710264220470329 # Parent 682f31b20894d2f4c6fef5650f13d7e385b47ddb [project @ 1996-02-12 03:07:39 by jwe] diff -r 682f31b20894 -r 333923894848 liboctave/CmplxAEPBAL.cc --- a/liboctave/CmplxAEPBAL.cc Mon Feb 12 02:26:07 1996 +0000 +++ b/liboctave/CmplxAEPBAL.cc Mon Feb 12 03:07:39 1996 +0000 @@ -50,44 +50,47 @@ int ComplexAEPBALANCE::init (const ComplexMatrix& a, const string& balance_job) { - int a_nc = a.cols (); + int n = a.cols (); - if (a.rows () != a_nc) + if (a.rows () != n) { (*current_liboctave_error_handler) ("AEPBALANCE requires square matrix"); return -1; } - int n = a_nc; - - // Parameters for balance call. - int info; int ilo; int ihi; - double *scale = new double [n]; - // Copy matrix into local structure. + Array scale (n); + double *pscale = scale.fortran_vec (); balanced_mat = a; - - char bal_job = balance_job[0]; + Complex *p_balanced_mat = balanced_mat.fortran_vec (); - F77_FCN (zgebal, ZGEBAL) (&bal_job, n, - balanced_mat.fortran_vec (), n, ilo, ihi, - scale, info, 1L, 1L); + char job = balance_job[0]; - // Initialize balancing matrix to identity. + F77_XFCN (zgebal, ZGEBAL, (&job, n, p_balanced_mat, n, ilo, ihi, + pscale, info, 1L, 1L)); - balancing_mat = Matrix (n, n, 0.0); - for (int i = 0; i < n; i++) - balancing_mat (i, i) = 1.0; + if (f77_exception_encountered) + (*current_liboctave_error_handler) ("unrecoverable error in zgebal"); + else + { + balancing_mat = Matrix (n, n, 0.0); + for (int i = 0; i < n; i++) + balancing_mat.elem (i, i) = 1.0; - F77_FCN (zgebak, ZGEBAK) (&bal_job, "R", n, ilo, ihi, scale, n, - balancing_mat.fortran_vec (), n, info, 1L, - 1L); + Complex *p_balancing_mat = balancing_mat.fortran_vec (); + + char side = 'R'; - delete [] scale; + F77_XFCN (zgebak, ZGEBAK, (&job, &side, n, ilo, ihi, pscale, n, + p_balancing_mat, n, info, 1L, 1L)); + + if (f77_exception_encountered) + (*current_liboctave_error_handler) ("unrecoverable error in zgebak"); + } return info; } diff -r 682f31b20894 -r 333923894848 liboctave/CmplxAEPBAL.h --- a/liboctave/CmplxAEPBAL.h Mon Feb 12 02:26:07 1996 +0000 +++ b/liboctave/CmplxAEPBAL.h Mon Feb 12 03:07:39 1996 +0000 @@ -59,6 +59,8 @@ return *this; } + ~ComplexAEPBALANCE (void) { } + ComplexMatrix balanced_matrix (void) const { return balanced_mat; } ComplexMatrix balancing_matrix (void) const { return balancing_mat; } diff -r 682f31b20894 -r 333923894848 liboctave/dbleAEPBAL.cc --- a/liboctave/dbleAEPBAL.cc Mon Feb 12 02:26:07 1996 +0000 +++ b/liboctave/dbleAEPBAL.cc Mon Feb 12 03:07:39 1996 +0000 @@ -49,44 +49,47 @@ int AEPBALANCE::init (const Matrix& a, const string& balance_job) { - int a_nc = a.cols (); + int n = a.cols (); - if (a.rows () != a_nc) + if (a.rows () != n) { (*current_liboctave_error_handler) ("AEPBALANCE requires square matrix"); return -1; } - int n = a_nc; - - // Parameters for balance call. - int info; int ilo; int ihi; - double *scale = new double [n]; - // Copy matrix into local structure. + Array scale (n); + double *pscale = scale.fortran_vec (); balanced_mat = a; - - char bal_job = balance_job[0]; + double *p_balanced_mat = balanced_mat.fortran_vec (); - F77_FCN (dgebal, DGEBAL) (&bal_job, n, - balanced_mat.fortran_vec (), n, ilo, ihi, - scale, info, 1L, 1L); + char job = balance_job[0]; - // Initialize balancing matrix to identity. + F77_XFCN (dgebal, DGEBAL, (&job, n, p_balanced_mat, n, ilo, ihi, + pscale, info, 1L, 1L)); - balancing_mat = Matrix (n, n, 0.0); - for (int i = 0; i < n; i++) - balancing_mat.elem (i ,i) = 1.0; + if (f77_exception_encountered) + (*current_liboctave_error_handler) ("unrecoverable error in dgebal"); + else + { + balancing_mat = Matrix (n, n, 0.0); + for (int i = 0; i < n; i++) + balancing_mat.elem (i ,i) = 1.0; - F77_FCN (dgebak, DGEBAK) (&bal_job, "R", n, ilo, ihi, scale, n, - balancing_mat.fortran_vec (), n, info, 1L, - 1L); + double *p_balancing_mat = balancing_mat.fortran_vec (); + + char side = 'R'; - delete [] scale; + F77_XFCN (dgebak, DGEBAK, (&job, &side, n, ilo, ihi, pscale, n, + p_balancing_mat, n, info, 1L, 1L)); + + if (f77_exception_encountered) + (*current_liboctave_error_handler) ("unrecoverable error in dgebak"); + } return info; } diff -r 682f31b20894 -r 333923894848 liboctave/dbleAEPBAL.h --- a/liboctave/dbleAEPBAL.h Mon Feb 12 02:26:07 1996 +0000 +++ b/liboctave/dbleAEPBAL.h Mon Feb 12 03:07:39 1996 +0000 @@ -59,6 +59,8 @@ return *this; } + ~AEPBALANCE (void) { } + Matrix balanced_matrix (void) const { return balanced_mat; } Matrix balancing_matrix (void) const { return balancing_mat; } diff -r 682f31b20894 -r 333923894848 liboctave/dbleGEPBAL.cc --- a/liboctave/dbleGEPBAL.cc Mon Feb 12 02:26:07 1996 +0000 +++ b/liboctave/dbleGEPBAL.cc Mon Feb 12 03:07:39 1996 +0000 @@ -73,14 +73,15 @@ int n = a_nc; - // Parameters for balance call. - int info; int ilo; int ihi; - double *cscale = new double [n]; - double *cperm = new double [n]; + + Array cscale (n); + double *pcscale = cscale.fortran_vec (); + Matrix wk (n, 6, 0.0); + double *pwk = wk.fortran_vec (); // Back out the permutations: // @@ -102,23 +103,18 @@ balanced_a_mat = a; balanced_b_mat = b; - // Initialize balancing matrices to identity. - - left_balancing_mat = Matrix (n, n, 0.0); - for (int i = 0; i < n; i++) - left_balancing_mat (i, i) = 1.0; - - right_balancing_mat = left_balancing_mat; + double *p_balanced_a_mat = balanced_a_mat.fortran_vec (); + double *p_balanced_b_mat = balanced_b_mat.fortran_vec (); // Check for permutation option. - char bal_job = balance_job[0]; + char job = balance_job[0]; - if (bal_job == 'P' || bal_job == 'B') + if (job == 'P' || job == 'B') { - F77_FCN (reduce, REDUCE) (n, n, balanced_a_mat.fortran_vec (), - n, balanced_b_mat.fortran_vec (), ilo, - ihi, cscale, wk.fortran_vec ()); + F77_XFCN (reduce, REDUCE, (n, n, p_balanced_a_mat, n, + p_balanced_b_mat, ilo, ihi, + pcscale, pwk)); } else { @@ -128,66 +124,102 @@ ihi = n; } - // Check for scaling option. - - if ((bal_job == 'S' || bal_job == 'B') && ilo != ihi) - { - F77_FCN (scaleg, SCALEG) (n, n, balanced_a_mat.fortran_vec (), - n, balanced_b_mat.fortran_vec (), ilo, - ihi, cscale, cperm, wk.fortran_vec ()); - } + if (f77_exception_encountered) + (*current_liboctave_error_handler) ("unrecoverable error in reduce"); else { - // Set scaling data to 0's. + Array cperm (n); + double *pcperm = cperm.fortran_vec (); + + // Check for scaling option. + + if ((job == 'S' || job == 'B') && ilo != ihi) + { + F77_XFCN (scaleg, SCALEG, (n, n, p_balanced_a_mat, n, + p_balanced_b_mat, ilo, ihi, + pcscale, pcperm, pwk)); + } + else + { + // Set scaling data to 0's. + + for (int i = ilo-1; i < ihi; i++) + { + cscale.elem (i) = 0.0; + wk.elem (i, 0) = 0.0; + } + } + + if (f77_exception_encountered) + (*current_liboctave_error_handler) ("unrecoverable error in scaleg"); + else + { + // Scaleg returns exponents, not values, so... + + for (int i = ilo-1; i < ihi; i++) + { + cscale.elem (i) = pow (2.0, cscale.elem (i)); + wk.elem (i, 0) = pow (2.0, -wk.elem (i, 0)); + } + + // Initialize balancing matrices to identity. + + left_balancing_mat = Matrix (n, n, 0.0); + for (int i = 0; i < n; i++) + left_balancing_mat (i, i) = 1.0; + + right_balancing_mat = left_balancing_mat; + + double *p_right_balancing_mat = right_balancing_mat.fortran_vec (); + double *p_left_balancing_mat = left_balancing_mat.fortran_vec (); - for (int tmp = ilo-1; tmp < ihi; tmp++) - { - cscale[tmp] = 0.0; - wk.elem (tmp, 0) = 0.0; + // Column permutations/scaling. + + char side = 'R'; + + F77_XFCN (dgebak, DGEBAK, (&job, &side, n, ilo, ihi, pcscale, + n, p_right_balancing_mat, n, info, + 1L, 1L)); + + if (f77_exception_encountered) + (*current_liboctave_error_handler) + ("unrecoverable error in dgebak"); + else + { + // Row permutations/scaling. + + side = 'L'; + + F77_XFCN (dgebak, DGEBAK, (&job, &side, n, ilo, ihi, pwk, + n, p_left_balancing_mat, n, + info, 1L, 1L)); + +#if 0 + // XXX FIXME XXX --- these four lines need to be added and + // debugged. GEPBALANCE::init will work without them, though, so + // here they are. + + if ((job == 'P' || job == 'B') && ilo != ihi) + { + F77_XFCN (gradeq, GRADEQ, (n, n, p_balanced_a_mat, n, + p_balanced_b_mat, ilo, ihi, + pcperm, pwk)); + } +#endif + + if (f77_exception_encountered) + (*current_liboctave_error_handler) + ("unrecoverable error in dgebak"); + else + { + // Transpose for aa = cc*a*dd convention... + + left_balancing_mat = left_balancing_mat.transpose (); + } + } } } - // Scaleg returns exponents, not values, so... - - for (int tmp = ilo-1; tmp < ihi; tmp++) - { - cscale[tmp] = pow (2.0, cscale[tmp]); - wk.elem (tmp, 0) = pow (2.0, -wk.elem (tmp, 0)); - } - - // Column permutations/scaling. - - F77_FCN (dgebak, DGEBAK) (&bal_job, "R", n, ilo, ihi, cscale, n, - right_balancing_mat.fortran_vec (), n, - info, 1L, 1L); - - // Row permutations/scaling. - - F77_FCN (dgebak, DGEBAK) (&bal_job, "L", n, ilo, ihi, - wk.fortran_vec (), n, - left_balancing_mat.fortran_vec (), n, - info, 1L, 1L); - - // XXX FIXME XXX --- these four lines need to be added and - // debugged. GEPBALANCE::init will work without them, though, so - // here they are. - -#if 0 - if ((bal_job == 'P' || bal_job == 'B') && ilo != ihi) - { - F77_FCN (gradeq, GRADEQ) (n, n, balanced_a_mat.fortran_vec (), - n, balanced_b_mat.fortran_vec (), ilo, - ihi, cperm, wk.fortran_vec ()); - } -#endif - - // Transpose for aa = cc*a*dd convention... - - left_balancing_mat = left_balancing_mat.transpose (); - - delete [] cscale; - delete [] cperm; - return info; } diff -r 682f31b20894 -r 333923894848 liboctave/dbleGEPBAL.h --- a/liboctave/dbleGEPBAL.h Mon Feb 12 02:26:07 1996 +0000 +++ b/liboctave/dbleGEPBAL.h Mon Feb 12 03:07:39 1996 +0000 @@ -66,6 +66,8 @@ return *this; } + ~GEPBALANCE (void) { } + Matrix balanced_a_matrix (void) const { return balanced_a_mat; } Matrix balanced_b_matrix (void) const { return balanced_b_mat; }