[Trilinos-Users] condition number of an Epetra_FECrsMatrix
Nate Roberts
nate at nateroberts.com
Wed Apr 17 09:52:25 MDT 2013
Pavel,
You are exactly correct: I have verified that MATLAB's cond() gives me 2640.2 for my example matrix: it's at least the same order of magnitude as what AztecOO was giving me. I'll investigate which norm is the more appropriate for my purposes, and have a look at Anasazi if I decide that I need something more accurate or in a different norm than what AztecOO is giving me.
Thanks for your help!
Nate
On Apr 17, 2013, at 10:31 AM, Pavel Jiránek <pavel.jiranek at gmail.com> wrote:
> Hi Nate,
>
> I believe that this is due to the fact that both the Matlab functions "condest" and "rcond" estimate the condition number of a given matrix with respect to the 1-norm while the Aztec's estimators (using an eigen- or singular-value analysis of the coefficient matrix from the Lanczos or Arnoldi process) estimate (most likely) the 2-norm. Condition numbers with respect to these norms are equivalent in the sense that one is a scalar multiple of the other with the scalar lying somewhere in the interval [1/n, n], where n is the matrix dimension (for square matrices).
>
> In order to verify your results from Aztec, you could try for a small matrix to use instead the function "cond" on the matrix converted to the full matrix (cond(full(A))) which computes the 2-norm condition number using the singular value decomposition. Also note that the quality of the estimate obtained from the coefficient matrix in Lanczos and Arnoldi process (used by Aztec) is strongly dependent on the "quality" of the projection (which depends on the convergence of the extreme eigen/singular values). You could also use probably some more "advanced" methods for eigenvalue problems in Anasazi (probably LOBPCG?) if you need an accurate result.
>
> Best wishes,
>
> Pavel
>
>
>
> 2013/4/17 Nate Roberts <nate at nateroberts.com>
> Hello all,
>
> I'm interested in studying the condition number of our global system matrix, which we store as an Epetra_FECrsMatrix. My approach thus far doesn't give accurate estimates, and I'd like advice on getting better ones.
>
> Based on some sample code I found in the AztecOO documentation, I wrote the following method:
>
> double Solution::conditionNumberEstimate( Epetra_LinearProblem & problem ) {
> AztecOOConditionNumber conditionEstimator;
> conditionEstimator.initialize(*problem.GetOperator());
>
> int maxIters = 40000;
> double tol = 1e-10;
> int status = conditionEstimator.computeConditionNumber(maxIters, tol);
> if (status!=0)
> cout << "status result from computeConditionNumber(): " << status << endl;
> double condest = conditionEstimator.getConditionNumber();
>
> return condest;
> }
>
> (The input Epetra_LinearProblem is initialized with the global system matrix.)
>
> For some smaller problems, I also tried exporting the matrix to file, and loading it into MATLAB, and using rcond() and/or condest() there. What I found was that the method above significantly underreports the condition number. For example, for a fairly dense 243x243 matrix (corresponding to a 2x2 finite element mesh), the method above reports 2611.05 as the condition number. MATLAB's condest and rcond both report 16716. (Decreasing the tolerance to 1e-12 and increasing maxIters to 100000, I instead get 2267.47 from the method above.)
>
> (It may be worth noting that the Epetra_FECrsMatrix ends up storing quite a few zeros. Of the 19237 entries stored, 5024 of them are zeros.)
>
> I'd really like to be able to get high-fidelity condition numbers for this matrix without going to MATLAB. I don't mind if it takes some computational effort—I'm not using the condition number during my real solves, but studying the effect of various parameters (e.g. choice of basis functions) on conditioning.
>
> Any advice would be greatly appreciated!
>
> Thanks,
> Nate
>
> _______________________________________________
> Trilinos-Users mailing list
> Trilinos-Users at software.sandia.gov
> http://software.sandia.gov/mailman/listinfo/trilinos-users
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: https://software.sandia.gov/pipermail/trilinos-users/attachments/20130417/15e30a82/attachment.html
More information about the Trilinos-Users
mailing list