[Trilinos-Users] Inverse of Serial Dense Matrix
Nico Schlömer
nico.schloemer at gmail.com
Fri Feb 17 08:21:36 MST 2012
Hi Paul,
there are specific classes for it; check out
http://trilinos.sandia.gov/packages/docs/dev/packages/teuchos/doc/html/classTeuchos_1_1SerialDenseSolver.html
http://trilinos.sandia.gov/packages/docs/dev/packages/teuchos/doc/html/classTeuchos_1_1SerialSpdDenseSolver.html
Here's how I used SerialSpdDenseSolver:
Quick & dirty:
================= *snip* =================
Teuchos::SerialSpdDenseSolver<int,double> solver;
TEUCHOS_ASSERT_EQUALITY( 0, solver.setMatrix( A ) );
TEUCHOS_ASSERT_EQUALITY( 0, solver.setVectors( alpha, rhs ) );
TEUCHOS_ASSERT_EQUALITY( 0, solver.solve() );
================= *snap* =================
Advanced:
================= *snip* =================
Teuchos::SerialSpdDenseSolver<int,double> solver;
TEUCHOS_ASSERT_EQUALITY( 0, solver.setMatrix( A ) );
TEUCHOS_ASSERT_EQUALITY( 0, solver.setVectors( alpha, rhs ) );
if ( solver.shouldEquilibrate() )
{
TEUCHOS_ASSERT_EQUALITY( 0, solver.equilibrateRHS() );
#if TRILINOS_MAJOR_MINOR_VERSION > 100803
TEUCHOS_ASSERT_EQUALITY( 0, solver.equilibrateMatrix() );
TEUCHOS_ASSERT_EQUALITY( 0, solver.solve() );
TEUCHOS_ASSERT_EQUALITY( 0, solver.unequilibrateLHS() );
#else
// A bug in Trilinos<=10.8.3 makes unequilibrateLHS() fail.
// Work around by doing the unequilibration manually.
// Note that this relies on the scaling being
// s(i) = 1 / sqrt(A(i,i)).
Teuchos::RCP<DoubleVector> diagA =
Teuchos::rcp( new DoubleVector( numEdges ) );
for ( int k=0; k<numEdges; k++ )
(*diagA)(k) = (*A)(k,k);
TEUCHOS_ASSERT_EQUALITY( 0, solver.equilibrateMatrix() );
TEUCHOS_ASSERT_EQUALITY( 0, solver.solve() );
for ( int k=0; k<numEdges; k++ )
(*alpha)(k) = (*alpha)(k) / sqrt( (*diagA)(k) );
#endif
}
else
{
TEUCHOS_ASSERT_EQUALITY( 0, solver.solve() );
}
================= *snap* =================
If you're running on a recent Trilinos version, you can forget about
the #else part of course.
Cheers,
Nico
On Fri, Feb 17, 2012 at 2:54 PM, Paul Dionne <pjd at cfdrc.com> wrote:
> What is the best method to obtain the inverse of a serial dense matrix?
> I think I'm misunderstanding the ApplyInverse() function. I also tried
> setting up a SerialDenseSolver problem and calling that invert()
> function, but that left the matrix unchanged.
>
> Paul Dionne
>
> _______________________________________________
> Trilinos-Users mailing list
> Trilinos-Users at software.sandia.gov
> http://software.sandia.gov/mailman/listinfo/trilinos-users
>
More information about the Trilinos-Users
mailing list