[Trilinos-Users] Beginner's problems with PyTrilinos and Anasazi

Marko Budisic mbudisic at engineering.ucsb.edu
Wed Dec 17 12:25:25 MST 2008


Hi Heidi,

yes, my test problem works fine now. Thank you and Chris for your help!
I am glad that my problem has done some greater good too :)

If I encounter any other problems along the way, I'll post again, but
for now I think I am set.

Best,
Marko

On Wed, 2008-12-17 at 12:15 -0700, Thornquist, Heidi K wrote:
> Hi Marko,
> 
> Thanks, I think you found a bug in one of our orthogonalization
> methods :)
> The initial vector that you pass in is initialized to zero, which is
> normally caught by 
> our orthogonalization methods, so can you generate the init vector by:
> 
> blocksize = 1
> storage = Epetra.MultiVector( W.RowMap() , blocksize)
> storage.Random();
> eigs.setInitVec(storage)
> 
> Let me know if this helps.
> 
> Thanks,
> Heidi
> 
> 
> On 12/16/08 11:28 AM, "Marko Budisic" <mbudisic at engineering.ucsb.edu>
> wrote:
> 
>         Hi Heidi,
>         
>         thank you for your reply.
>         
>         Indeed, I have intended to create a 5x5 matrix and the code
>         snippet had
>         a faulty code, probably left over from when I was trying to
>         play with
>         parameters to see if I understood them properly. The line
>         
>         numglobals = nx*nx
>         
>         should be changed to
>         
>         numglobals = nx
>         
>         to get a matrix with 5 rows. The matrix population code is
>         complete.
>         However, this does not resolve the original problem. For
>         additional
>         information, here are eigenvalues of the matrix that I would
>         like to
>         compute (as computed by numpy and numpy.linalg):
>         
>         1.00000000+0.j
>         -0.01363536+0.59362319j
>         -0.01363536-0.59362319j
>         0.06893738+0.j
>         0.33333333+0.j
>         
>         Thank you for your help.
>         
>         Best,
>         Marko
>         
>         
>         On Tue, 2008-12-16 at 11:09 -0700, Thornquist, Heidi K wrote:
>         > Hi Marko,
>         >
>         > Thanks for sending a copy of your code.  Upon inspection of
>         the
>         > attached code snippet, it looks like you
>         > are creating a matrix with 25 rows and columns and then only
>         filling
>         > the primary 5x5 block.  This gives
>         > the matrix a nullspace of dimension 20.  Is this the matrix
>         you meant
>         > to generate, or did you omit some
>         > of the matrix population for brevity?  If so, could you
>         include all
>         > the matrix population code.
>         >
>         > Thanks,
>         > Heidi
>         >
>         > ---
>         > Heidi K. Thornquist
>         > Senior Member of Technical Staff
>         > Electrical & Microsystem Modeling
>         > Sandia National Laboratories
>         > P.O. Box 5800, MS 0316
>         > Albuquerque, NM  87185-0316
>         > Phone: (505)284-8426
>         > Fax:   (505)284-5451
>         >
>         >
>         >
>         >
>         > On 12/15/08 9:47 PM, "Marko Budisic"
>         <mbudisic at engineering.ucsb.edu>
>         > wrote:
>         >
>         >         Dear Trilinos list,
>         >        
>         >         I would like to use PyTrilinos with Anasazi to
>         compute
>         >         eigenvectors of
>         >         large non-symmetric sparse matrices. In order to
>         familiarize
>         >         myself with
>         >         both packages, I tried to build a small example for
>         computing
>         >         1-eigenvector of a probability transition matrix for
>         a random
>         >         walk on a
>         >         graph. The solver I tried to use was
>         BlockKrylovSchur.
>         >        
>         >         I managed to successfuly build the matrix in
>         question and even
>         >         perform
>         >         some basic operations on it (I checked by printing
>         the
>         >         matrix). I also
>         >         managed to pass the setProblem() sanity test with
>         the
>         >         Eigenproblem
>         >         object. However, when I try to invoke the solver,
>         the GEES
>         >         subroutine
>         >         fails which causes the entire code to fail.
>         >        
>         >         This is the traceback from Python:
>         >         ******************
>         >         Traceback (most recent call last):
>         >           File "src/anasazitest.py", line 56, in <module>
>         >             info = mngr.solve()
>         >           File
>         >
>                 "/usr/lib/python2.5/site-packages/PyTrilinos/Anasazi.py", line
>         >         6456, in solve
>         >             return
>         _Anasazi.BlockKrylovSchurSolMgrEpetra_solve(*args)
>         >
>                 SystemError: /home/mbudisic/downloads/trilinos-9.0.1/packages/anasazi/src/AnasaziBlockKrylovSchur.hpp:1432:
>         >        
>         >         Throw number = 1
>         >        
>         >         Throw test that evaluated to true: info != 0
>         >        
>         >         Anasazi::BlockKrylovSchur::computeSchurForm(): GEES
>         returned
>         >         info != 0.
>         >         *******************
>         >        
>         >         I am guessing that I am mismatching some options,
>         possibly
>         >         related to
>         >         block size and sizes of auxiliary vectors but I
>         don't know how
>         >         to
>         >         correct the problem. Below is the relevant source
>         code.
>         >        
>         >         If anyone could direct me in correcting this
>         problem, I would
>         >         be very
>         >         thankful. Also, if you have any other suggestions,
>         based on my
>         >         code,
>         >         those would also be very much appreciated.
>         >        
>         >         Best,
>         >         Marko Budisic
>         >        
>         >         *********************************
>         >         # epetra
>         >         nx = 5
>         >         numglobals = nx*nx
>         >         mycomm = Epetra.SerialComm()
>         >         mymap = Epetra.Map( numglobals, 0, mycomm )
>         >         W = Epetra.CrsMatrix(Epetra.Copy, mymap, 4)
>         >        
>         >         # populate matrix
>         >         W.InsertGlobalValues( 0, [1.,2.,3.,2.], [0,1,2,4] )
>         >         W.InsertGlobalValues( 1, [1.,2.], [1,2] )
>         >         W.InsertGlobalValues( 2, [1.,2.], [2,3] )
>         >         W.InsertGlobalValues( 3, [3.,1.], [0,3] )
>         >         W.InsertGlobalValues( 4, [2.,1.], [3,4] )
>         >         W.FillComplete()
>         >        
>         >         # normalize the matrix
>         >         d = Epetra.Vector( W.RowMap() )
>         >         W.InvRowSums(d)
>         >         W.LeftScale(d)
>         >        
>         >         print W
>         >        
>         >         eigs = Anasazi.BasicEigenproblem()
>         >         eigs.setA( W )
>         >         eigs.setHermitian(False)
>         >        
>         >         Nev = 1
>         >         eigs.setNEV( Nev )
>         >        
>         >         blocksize = 1
>         >         storage = Epetra.MultiVector( W.RowMap() ,
>         blocksize)
>         >         eigs.setInitVec(storage)
>         >        
>         >         print "Problem set:", eigs.setProblem()
>         >        
>         >         myPL = Teuchos.ParameterList()
>         >         myPL.set("Which", "LM")
>         >         myPL.set("Block Size",blocksize)
>         >         myPL.set("Verbosity", Anasazi.Errors |
>         Anasazi.Warnings | \
>         >                      Anasazi.FinalSummary )
>         >        
>         >         mngr = Anasazi.BlockKrylovSchurSolMgr( eigs, myPL )
>         >         info = mngr.solve()
>         >         print "Problem solved:", info
>         >        
>         >         x = eigs.getSolution()
>         >         print x.Evecs()
>         >        
>         >         *********************************
>         >        
>         >        
>         >         _______________________________________________
>         >         Trilinos-Users mailing list
>         >         Trilinos-Users at software.sandia.gov
>         >
>                 http://software.sandia.gov/mailman/listinfo/trilinos-users
>         >        
>         >
>         > _______________________________________________
>         > Trilinos-Users mailing list
>         > Trilinos-Users at software.sandia.gov
>         > http://software.sandia.gov/mailman/listinfo/trilinos-users
>         
>         
>         _______________________________________________
>         Trilinos-Users mailing list
>         Trilinos-Users at software.sandia.gov
>         http://software.sandia.gov/mailman/listinfo/trilinos-users
>         
> 
> _______________________________________________
> 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