[Trilinos-Users] Anasazi report error

maireni at cimne.upc.edu maireni at cimne.upc.edu
Mon Dec 15 16:02:45 MST 2014


Sorry Alice,
let me be more clear

What is the size of your matrix? 27100,be
What is the "Block Size" you're using? 3
What is "Num Blocks"? 3
You think is not good idea to mixed num Blocks and number of Restart 
Blocks?
Regards

Nelson





El 2014-12-16 02:03, Alicia Klinvex escribió:
> I can't tell for sure from the code you provided, but it looks like
> you may have mixed up numBlocks and numRestartBlocks.
> 
> How big is your matrix?
> 
> On Mon, Dec 15, 2014 at 9:56 AM, <maireni at cimne.upc.edu> wrote:
> 
>> Hi Alice.
>> Thanks for answer.
>> Here is the data
>> 
>> double tol      = 1.0e-8; // convergence tolerance
>> int nev         = eig_p->number_eig; // number of eigenvalues
>> for which to solve
>> int blockSize   = 3; // block size (number of eigenvectors
>> processed at once)
>> int numBlocks   = 3 * nev / blockSize; // restart length
>> int maxRestarts = 5; // maximum number of restart cycles
>> int numRestartBlocks = 8*nev/blockSize;
>> std::string whenToShift = "Always";
>> 
>> Thanks again.
>> 
>> Nelson
>> 
>> El 2014-12-15 22:42, Alicia Klinvex escribió:
>> 
>> Hello,
>> 
>> That error is trying to tell you that there are too many vectors in
>> your subspace, and block Krylov Schur can't orthonormalize them. 
>> Example: your matrix is 8x8, block size is 2, and the maximum
>> number
>> of blocks is 5.  That's not going to work, because BKS is going to
>> try to orthonormalize 10 vectors (2*5), but the vectors are of
>> length
>> 8, so it will fail.
>> 
>> Three things:
>> * What is the size of your matrix?
>> * What is the "Block Size" you're using?
>> * What is "Num Blocks"?
>> 
>> - Alicia
>> 
>> On Mon, Dec 15, 2014 at 4:18 AM, <maireni at cimne.upc.edu> wrote:
>> 
>> Hi Alicia,
>> I have some good news for using Anasazi solver to solve my eigen
>> problems.
>> As you know, i am using the AMLS method, and the key here is
>> that a system (size of matrix), let say a half million dofs, just i
>> can reduced to some thousand
>> dofs.
>> 
>> However, in my final eigen problem my new mass matrix is really
>> populated,
>> for taking a number, the dimension is 27100x27100 and just 43% are
>> non zeros.
>> My new stiff is diagonal matrix. When i call Anasazi with Amesos
>> (using taucs solver)
>> i have the following report:
>> 
>> terminate called after throwing an instance of
>> 'std::invalid_argument'
>>   what(): 
>> 
>> 
> /usr/local/trilinos-11.12.1/include/AnasaziBlockKrylovSchurSolMgr.hpp:338:
>> 
>> Throw number = 1
>> 
>> Throw test that evaluated to true:
>> static_cast<ptrdiff_t>(_numBlocks)*_blockSize >
>> MVText::GetGlobalLength(*_problem->getInitVec())
>> 
>> Anasazi::BlockKrylovSchurSolMgr: Potentially impossible
>> orthogonality requests. Reduce basis size.
>> 
>> How i can make it Anasazi works with the characteristic of matrix
>> say before.
>> (K->Diagonal and M->sparse matrix but do not so sparse or do not so
>> dense).
>> Thanks a lot for your help.
>> 
>> Regards
>> 
>> Nelson
>> 
>> El 2014-11-18 01:35, Alicia Klinvex escribió:
>> Hello,
>> 
>> I'm an external collaborator, so I really don't know much about
>> AztecOO.  (The only things I wrote in Trilinos are TraceMin and
>> some
>> of the things TraceMin depends on.)  If you have questions about
>> TraceMin, I would be happy to answer them.  Otherwise, you're more
>> likely to get a useful response by emailing the trilinos-users
>> list. 
>> When you email the trilinos-users list, your message will reach ALL
>> of
>> the developers (including me).
>> 
>> Hopefully one of the AztecOO developers will be able to answer your
>> question :-)
>> 
>> Best wishes,
>> Alicia
>> 
>> ---------- Forwarded message ----------
>> From: <maireni at cimne.upc.edu>
>> Date: Mon, Nov 17, 2014 at 11:11 AM
>> Subject: Re: Krylov-Schur method with AztecOO
>> To: Alicia Klinvex <aklinvex at purdue.edu>
>> 
>> Hi Alicia,
>> 
>>  I understand how work Epetra matrix and i was able to
>>  calculate eigen values with  block Krylov-Schur method
>>  using Amesos. However, changing to AztecOO i have a "nan"
>>  as you see above. I dont know, but in this case
>>  K and M are sparse stiffness and mass matrix respectively.
>>  M is diagonal and have some zero terms in diagonal (rotational
>> Dofs).
>> 
>>  I don do  understand what is the problem  exactly.
>>  Do you you have any idea about this mistake?.
>>  Or just using AztecOO is limited to some special matrices,
>>  due is a iterative solver.
>> 
>>  Regards and best wishes
>> 
>>  Nelson
>> 
>>                      iter:    0         
>>  residual =
>> 1.000000e+00
>>                      iter:    1         
>>  residual =
>> -nan
>> 
>>         
>> ***************************************************************
>> 
>>          Warning: a breakdown in this method
>>          has occurred and solution has not converged.
>> 
>>          Solver:                 bicgstab
>>          number of iterations:   1
>> 
>>          Recursive residual =        -nan
>> 
>>  Anasazi::BlockKrylovSchurSolMgr::solve() caught unexpected
>> exception
>> from Anasazi::BlockKrylovSchur::iterate() at iteration 1
>>  /usr/local/trilinos-11.12.1/include/AnasaziEpetraAdapter.hpp:1401:
>> 
>>  Throw number = 1
>> 
>>  Throw test that evaluated to true: ret != 0
>> 
>> 
>  Anasazi::OperatorTraits<double,Epetra_Multivector,Epetra_Operator>::Apply():
>> Error in Epetra_Operator::Apply(). Code -2
>>  Anasazi::BlockKrylovSchurSolMgr::solve() returning Unconverged
>> with
>> no solutions.
>>  Anasazi::EigensolverMgr::solve() returned unconverged.
>> 
>>  El 2014-11-13 02:49, Alicia Klinvex escribió:
>> 
>> The convergence rate of each eigenpair is based on the quantity
>> lambda_i / lambda_{s+1}, where s is the subspace dimension you
>> selected.  TraceMin by default selects 2*the number of desired
>> eigenpairs.  If the distribution of eigenvalues is bad, then it
>> can
>> take many iterations.  There are still things we can do to fix it,
>> like increasing the subspace dimension or using shifts.
>> 
>> Looking at that quantity lambda_i / lambda_{s+1} will tell you a
>> lot
>> about why TraceMin is behaving this way.
>> 
>> On Tue, Nov 11, 2014 at 10:50 AM, <maireni at cimne.upc.edu> wrote:
>> 
>> Hi,
>> I am using TraceMin for solve the problem.
>> My problem is symmetric and sparse.
>> However when i use lapack(is time consuming because is for dense
>> matrix)
>> FEAST or for example , TraceMin takes time for converge.
>> When i look for 4 eigen-values TraceMin is fast but when i
>> increase
>> the number of eigen-values
>> to 20 TraceMin make a lot iteration. I was thinking that maybe i
>> have a wrong parameter.
>> I dont know.
>> 
>> Regards
>> 
>> Thanks
>> 
>> Nel
>> 
>> El 2014-11-11 23:15, Alicia M Klinvex escribió:
>> 
>> Hello,
>> 
>> I'm a little confused by the fact that you're using LAPACK and
>> TraceMin on the same problem.  Is this a dense or sparse
>> problem?
>> 
>> Yes, all of the Anasazi solvers should work with Tpetra as well
>> as
>> Epetra.
>> 
>> Best wishes,
>> Alicia
>> 
>> ----- Original Message -----
>> From: maireni at cimne.upc.edu
>> To: "Alicia M Klinvex" <aklinvex at purdue.edu>
>> Sent: Tuesday, November 11, 2014 10:09:07 AM
>> Subject: Re: TraceMin-Davidson failurem
>> 
>> Hi Alice,
>> 
>> It is works now. I tested and give the similar values
>> of lapack eigen solver. However, why if  increase the number
>> of eigen-values to looking for (at instance from 4 to 20
>> eigen-values),
>> it is take more time.
>> My toler is 1E-6.
>> 
>> I will try other eigen-solver too. I suppose Tpetra
>> can be use in the others eigen-solver.
>> 
>> Regards and thanks again
>> 
>> Nel
>> 
>> El 2014-11-11 01:51, Alicia M Klinvex escribió:
>> Hello,
>> 
>> Thank you so much for giving our newest solver a try!
>> 
>> Anasazi calls a function called isFillComplete, which tells us
>> whether
>> you finished putting things in a matrix.  If you did NOT finish
>> putting things in the matrix, it breaks and outputs the error
>> message
>> you saw.
>> 
>> You may be thinking, "but I WAS finished modifying the entries of
>> the
>> matrix!"  You need to let Anasazi know that by calling the
>> Tpetra
>> function fillComplete
>> 
>> 
> (http://trilinos.org/docs/dev/packages/tpetra/doc/html/classTpetra_1_1CrsMatrix.html#a27c06cd71c7542742a731a6cc91180d2
>> [1]
>> [3]
>> 
>> [3]
>> 
>> [1]).
>>  In your code, this will look like K->fillComplete() and
>> M->fillComplete().
>> 
>> Does that answer your question?
>> 
>> Best wishes,
>> Alicia
>> 
>> ----- Original Message -----
>> From: maireni at cimne.upc.edu
>> To: "Alicia M Klinvex" <aklinvex at purdue.edu>
>> Sent: Monday, November 10, 2014 2:54:46 AM
>> Subject: Re: [Trilinos-Users] Returned mail: see transcript for
>> details
>> 
>> Hi Alicia,
>> At first so sorry for reply you a lot time for questions.
>> I am using the example you defined in the list that use
>> TraceMin-Davidson to solve a generalized eigenvalue problem
>> using Tpetra (and reading file).
>> 
>> The thing is i can not yet run my eigen problem.
>> In fact i have this mistakes and no reported in google
>> exist.  Here it is.
>> 
>> Teuchos::ArrayRCP<unsigned
> 
> 
> long>{ptr=0x46f9130,lowerOffset=0,upperOffset=8238,size=8239,node=0x4d6c7c0,strong_count=2,weak_count=0}
> 
>>>   >> Allocating X_
>>>   >> Allocating KX_
>>>   >> Allocating MX_
>>>   >> Allocating R_
>>>   >> Allocating V_
>>>   >> Allocating KV_
>>>   >> Allocating MV_
>>>   >> done allocating.
>>> Computing MV
>>> terminate called after throwing an instance of
>>> 'std::runtime_error'
>>>    what():
>>> /usr/local/trilinos-11.12.1/include/Tpetra_CrsMatrix_def.hpp:3517:
>>> 
>>> Throw number = 1
>>> 
>>> Throw test that evaluated to true: ! isFillComplete ()
>>> 
>>> What do you think i am doing wrong. I attached the code function
>>> if you have time to check it for a while. I really interesting
>>> about
>>> you tell me before, because i want to use any eigen-solver you
>>> defined
>>> in anasazi.
>>> 
>>> Thanks a lot
>>> 
>>> Nelson
>>> 
>>> El 2014-11-06 22:20, Alicia M Klinvex escribió:
>>> Hello,
>>> 
>>> a) The matrix should be the total matrix.  If your matrix comes
>>> from a
>>> matrix-market or harwell-boeing file, we do have a function that
>>> will
>>> read the symmetric matrix from a file and store it as general. 
>>> (I
>>> believe this is to optimize matrix-vector multiplications.)
>>> Alternatively, you can define your own matrix-vector
>>> multiplication
>>> operator, in which case Anasazi never needs to see your matrix.
>>> There's an example of this for TraceMin-Davidson.
>>> 
>>> b) Are you using BKS?  If you're not using block Krylov-Schur,
>>> you
>>> shouldn't need to create your own solver.  RTR and
>>> TraceMin/TraceMin-Davidson have solvers built-in.
>>> 
>>> c) I'm not sure I understand what you mean by the matrix being
>>> dense
>>> and sparse.  Also, are you looking for the smallest or largest
>>> eigenvalues?  Are you looking for interior or extreme
>>> eigenvalues?
>>> 
>>> By the way, one of the cool things about Anasazi is that all of
>>> the
>>> eigensolvers use the same interface.  If you want to try a
>>> couple
>>> of
>>> them, it only involves changing maybe 5 lines of code in your
>>> driver.
>>> 
>>> Best wishes,
>>> Alicia
>>> 
>>> PS: We should copy the user list in our correspondence, since
>>> some
>>> people are better equipped than I am to answer certain types of
>>> questions.  For instance, I am not an Amesos developer.  The
>>> Amesos
>>> developers can probably do a better job answering your questions
>>> about
>>> Amesos.
>>> 
>>> ----- Original Message -----
>>> From: maireni at cimne.upc.edu
>>> To: "Alicia M Klinvex" <aklinvex at purdue.edu>
>>> Sent: Thursday, November 6, 2014 7:05:34 AM
>>> Subject: Re: [Trilinos-Users] Returned mail: see transcript for
>>> details
>>> 
>>> Hello Alicia,
>>> 
>>> Just simple questions
>>> 
>>> a) When i use Anasazi with Tpetra in CSR format.
>>> Should be the total matrix even if it is symmetric or it is ok if
>>> i
>>> take
>>> for Anasazi the upper or lower one?. How to do that?
>>> 
>>> b) I see Amesos accept now Pardiso solver. I use normally this
>>> Solver.
>>> My matrices is sparse and lapack  is call when i use anasazi.
>>> How
>>> can
>>> i
>>> set amesos for use Pardiso or SuperLU(shared Memory). I see some
>>> macros
>>> but, i think should be another way no?
>>> 
>>> c) The problem i have is symmetric dense and sparse. According
>>> to your experience with is the fastest method in Anasazi?
>>> 
>>> Sorry again for ask you a lot, but i really excited with
>>> Trilinos.
>>> 
>>> Regards
>>> 
>>> Nelson
>>> 
>>> El 2014-11-06 10:38, Alicia M Klinvex escribió:
>>> Hello,
>>> 
>>> Have you tried contacting the SLEPC team about your problems?
>>> They're
>>> a helpful group as well.  (I'm obviously not trying to tell you
>>> not
>>> to
>>> use Trilinos; I just want you to be aware of the packages that
>>> can
>>> solve your specific problem.  We don't currently support finding
>>> all
>>> eigenvalues in an interval.  I should have that complete by the
>>> end
>>> of
>>> the month though.)
>>> 
>>> FEAST is slow because of its reliance on direct solvers such as
>>> PARDISO and all of the factorizations it requires.  You may want
>>> to
>>> contact the FEAST developers to see whether it would be
>>> appropriate
>>> to
>>> reduce the number of contour integration points, which can speed
>>> up
>>> the code.
>>> 
>>> Anasazi does have very powerful tools, and I'm glad you've
>>> decided
>>> to
>>> give them a try. :-)
>>> 
>>> I'm not familiar with AMLS.  Are these small matrices dense or
>>> sparse?
>>>  Regarding the linear systems, if they are sparse, Amesos and
>>> Belos
>>> are both really nice packages to use, depending on what your
>>> matrix
>>> is
>>> like.  Amesos contains direct solvers, and Belos is for
>>> iterative
>>> solvers such as CG.
>>> 
>>> Best wishes, and please keep us updated!
>>> Alicia
>>> 
>>> ----- Original Message -----
>>> From: maireni at cimne.upc.edu
>>> To: "Alicia M Klinvex" <aklinvex at purdue.edu>
>>> Sent: Wednesday, November 5, 2014 8:57:55 PM
>>> Subject: Re: [Trilinos-Users] Returned mail: see transcript for
>>> details
>>> 
>>> Hi Alicia M Klinvex,
>>> Thanks for replying.
>>> I was trying to use SLEPC's Krylov Schur but not successful yet.
>>> I know the  FEAST eigen-solver and is integrated in my code and
>>> works.
>>> I call directly from MKL libraries, but is slow than other
>>> eigen-solvers,
>>> at least block Lanzcos method like BLZPACK. Another problem is
>>> ARPACK++,
>>> that when i call the AReig functions not works correctly. It use
>>> a
>>> format
>>> of CSC. My matrices is symmetric and could be lower or upper
>>> matrix.
>>> In
>>> any case,
>>> the eigen-values is not ok. I saw Anasazi is inspired and have
>>> better
>>> tools than ARPACK. This is why i
>>> try to use Anasazi.
>>> 
>>> Respect to size of problem, the idea is to use AMLS method. AMLS
>>> method
>>> is a condensation method
>>> that can convert a huge matrix, millions of dofs, to just a
>>> little
>>> thousands. In need to solve
>>> a lot of eigen-values problems, but this eigen problems is so
>>> small
>>> that
>>> the corresponding one.
>>> Also i need to solve a lot equations with a lot rigth hand side.
>>> I
>>> was
>>> thinking in use
>>> Pardiso, but trilnos have good tools and i want to explot too.
>>> 
>>> Allow me try what you talk me and then i let you know.
>>> Thanks for your help. Thanks a lot.
>>> 
>>> Nel
>>> 
>>> El 2014-11-06 01:27, Alicia M Klinvex escribió:
>>> Hello Nelson,
>>> 
>>> What you want is this constructor:
> 
> 
> http://trilinos.org/docs/r11.12/packages/tpetra/doc/html/classTpetra_1_1Map.html#aeb556362352e2ae62804601e36672a6c
> [4]
>  [4]
> 
>> [4]
>> 
>>> [2]
>>> 
>>> Map( global_size_t numGlobalElements, GlobalOrdinal indexBase,
>>> const
>>> Teuchos::RCP< const Teuchos::Comm< int > > & comm, LocalGlobal lg
>>> =
>>> GloballyDistributed, const Teuchos::RCP< Node > & node =
>>> defaultArgNode< Node >() )
>>> 
>>> Here's what this constructor expects:
>>> * The number of rows of your matrix
>>> * The index base (which is generally 0)
>>> * The communicator (which we will come back to)
>>> * Whether the matrix is locally replicated or globally
>>> distributed
>>> (which we will come back to)
>>> * The node type, which defines the type of architecture you are
>>> using
>>> 
>>> You may be thinking at this point, "I have to provide a
>>> communicator?
>>> But I don't want to use MPI!"  You can define a serial Teuchos
>>> communicator which never communicates, as in this example:
> 
> 
> http://trilinos.org/docs/r11.12/packages/tpetra/doc/html/lesson01_no_mpi_8cpp-example.html
> [5]
>  [5]
> 
>> [5]
>> 
>>> [3]
>>> 
>>> For a sequential example, it doesn't matter whether you specify
>>> that
>>> the matrix is locally replicated or globally distributed.  We
>>> can
>>> come
>>> back to this when you're ready to use MPI.
>>> 
>>> Regarding finding all eigenvalues in an interval, you have
>>> several
>>> options:
>>> * You can compute the number of eigenvalues in an interval using
>>> inertia computations (meaning sparse factorizations will be
>>> involved),
>>> then call any eigensolver you want.
>>> * You can try Eric Polizzi's FEAST eigensolver (which requires an
>>> estimate of the number of eigenvalues in an interval).  You can
>>> download it here: http://www.feast-solver.org/ [2] [1] [1] [4] 
>>> It
>>> 
>>> can
>>> also be
>>> found
>>> in certain versions of the Intel Math Kernel Library.
>>> * You can try SLEPC's Krylov Schur, which is supposed to support
>>> this
>>> feature.  Just as a warning, I've personally never gotten it to
>>> work
>>> successfully.
>>> 
>>> We are currently working on adding that functionality (computing
>>> all
>>> eigenvalues in an interval) to Anasazi.  Just out of curiosity,
>>> roughly how many eigenvalues are there in your interval of
>>> interest?
>>> Also, how big is your matrix?
>>> 
>>> Best wishes,
>>> Alicia
>>> 
>>> ----- Original Message -----
>>> From: maireni at cimne.upc.edu
>>> To: "Alicia M Klinvex" <aklinvex at purdue.edu>
>>> Sent: Wednesday, November 5, 2014 1:08:00 AM
>>> Subject: Re: [Trilinos-Users] Returned mail: see transcript for
>>> details
>>> 
>>> Hi Alicia M Klinvex,
>>> 
>>> Thanks for reply.
>>> In my case, as in research area, i compile Trilinos for
>>> not running in MPI process, just in serial way. I know i can not
>>> exploits
>>> all tools of Trilinos, but is just for research. It have really
>>> interesting tools.
>>> When the code and the idea is working, migration to MPI is really
>>> well.
>>> 
>>> For the case of column map and row map i need your help.
>>> At the moment, just a class matrix works fine in Ananazi is ok
>>> for
>>> me.
>>> Epectra or Tpectra is ok. I am using now Tpectra and i see i can
>>> define
>>> a typedef
>>> like this
>>> 
>>> typedef double scalar_type;
>>> typedef Tpetra::CrsMatrix<scalar_type> TpetraCrsMatrix;
>>> 
>>> the other parameters is a default values, so i assumed with they
>>> Anasazi works. My dude is in how to contruct the maps for serial
>>> code.
>>> 
>>> Another thing i would like to ask you is about the eigen-solver.
>>> Anasazi can provide me the number of desired eigenvalues,
>>> but is possible it give to me an interval of eigen-values
>>> i.e [0, fecuency cuttof]?
>>> 
>>> I really appreciate your help.
>>> 
>>> Best wishes,
>>> 
>>> Nelson
>>> 
>>> El 2014-11-04 21:57, Alicia M Klinvex escribió:
>>> Hello Nelson,
>>> 
>>> My name is Alicia, and I am an Anasazi developer.  I'm glad
>>> you've
>>> decided to try Anasazi :-)
>>> 
>>> I assume what you're asking is how to take your three CSR arrays
>>> and
>>> turn them into a CrsMatrix that you can pass to an Anasazi
>>> eigensolver.  If you're using Tpetra, the constructor you're
>>> looking
>>> for is this:
> 
> 
> http://trilinos.org/docs/r11.10/packages/tpetra/doc/html/classTpetra_1_1CrsMatrix.html#a139a91ad9430f012bfb3bcb8fe413ae1
> [6]
>  [6]
> 
>> [6]
>> 
>>> [5]
>>> (If you're using Epetra instead, or want to know the difference
>>> between Epetra and Tpetra, let me know.)
>>> 
>>> Note that this constructor will expect the row map and column
>>> map,
>>> which define how the data is distributed over the MPI processes.
>>> If
>>> you don't understand how to create those maps, please let me know
>>> and
>>> I'll help you out.
>>> 
>>> Does that answer your question?
>>> 
>>> Best wishes,
>>> Alicia Klinvex
>>> 
>>> ----- Original Message -----
>>> From: maireni at cimne.upc.edu
>>> To: trilinos-users at software.sandia.gov
>>> Sent: Tuesday, November 4, 2014 3:41:49 AM
>>> Subject: Re: [Trilinos-Users] Returned mail: see transcript for
>>> details
>>> 
>>> Dear Heidi Thornquist,
>>> 
>>>   Allow me introduce myself first. My name is Nelson Lafontaine,
>>>   i am a PhD student from UPC, Barcelona, Spain. Now i am
>>> Beijing
>>> China
>>>   working with the Psor Chen Pu in the Automatic Multilevel
>>> Substructure
>>>   Methods (AMLS).
>>>   I see you are a Anasazi developer and i write to you because i
>>> am
>>>   interested in use
>>>   this library. My installation of Trilinos (trilinos-11.12.1)
>>> was
>>>   successfully,
>>>   even  i can run some test examples that you are define in
>>> Anasazi
>>> files.
>>> 
>>>   However, like you create your own class of matrices and
>>> vector,
>>>   my data i can provided to solve my eigen problem is in the
>>> form
>>> of
>>>   a Compressed Sparse Row (CSR). My question  is if exist a way
>>> to
>>> call
>>>   the eigen-solver defined in Anasazi just providing the data of
>>>   compressed
>>>   Sparse row format, or if is possible to create the matrices
>>> types
>>>   defined
>>>   in Anasazi just using the CRS format.
>>> 
>>>   Thanks in advance for your time.
>>> 
>>>   Nelson
>>> _______________________________________________
>>> Trilinos-Users mailing list
>>> Trilinos-Users at software.sandia.gov
>>> https://software.sandia.gov/mailman/listinfo/trilinos-users [3]
>>> [2]
>>> [2]
>>> [6]
>> 
>> Links:
>> ------
>> [1]
> 
> 
> http://trilinos.org/docs/dev/packages/tpetra/doc/html/classTpetra_1_1CrsMatrix.html#a27c06cd71c7542742a731a6cc91180d2
> [1]
>  [3]
> 
>> [3]
>> [2]
> 
> 
> http://trilinos.org/docs/r11.12/packages/tpetra/doc/html/classTpetra_1_1Map.html#aeb556362352e2ae62804601e36672a6c
> [4]
>  [4]
> 
>> [4]
>> [3]
> 
> 
> http://trilinos.org/docs/r11.12/packages/tpetra/doc/html/lesson01_no_mpi_8cpp-example.html
> [5]
>  [5]
> 
>> [5]
>> [4] http://www.feast-solver.org/ [2] [1] [1]
>> [5]
> 
> 
> http://trilinos.org/docs/r11.10/packages/tpetra/doc/html/classTpetra_1_1CrsMatrix.html#a139a91ad9430f012bfb3bcb8fe413ae1
> [6]
>  [6]
> 
>> [6]
>> [6] https://software.sandia.gov/mailman/listinfo/trilinos-users [3]
>> [2]
>> [2]
> 
>   Links:
>   ------
>   [1] http://www.feast-solver.org/ [2] [1]
>   [2] https://software.sandia.gov/mailman/listinfo/trilinos-users [3]
> [2]
>   [3]
> 
> 
> http://trilinos.org/docs/dev/packages/tpetra/doc/html/classTpetra_1_1CrsMatrix.html#a27c06cd71c7542742a731a6cc91180d2
> [1]
>  [3]
>   [4]
> 
> 
> http://trilinos.org/docs/r11.12/packages/tpetra/doc/html/classTpetra_1_1Map.html#aeb556362352e2ae62804601e36672a6c
> [4]
>  [4]
>   [5]
> 
> 
> http://trilinos.org/docs/r11.12/packages/tpetra/doc/html/lesson01_no_mpi_8cpp-example.html
> [5]
>  [5]
>   [6]
> 
> 
> http://trilinos.org/docs/r11.10/packages/tpetra/doc/html/classTpetra_1_1CrsMatrix.html#a139a91ad9430f012bfb3bcb8fe413ae1
> [6]
>  [6]
> 
>  Links:
>  ------
>  [1] http://www.feast-solver.org/ [2]
>  [2] https://software.sandia.gov/mailman/listinfo/trilinos-users [3]
>  [3]
> 
> http://trilinos.org/docs/dev/packages/tpetra/doc/html/classTpetra_1_1CrsMatrix.html#a27c06cd71c7542742a731a6cc91180d2
> [1]
>  [4]
> 
> http://trilinos.org/docs/r11.12/packages/tpetra/doc/html/classTpetra_1_1Map.html#aeb556362352e2ae62804601e36672a6c
> [4]
>  [5]
> 
> http://trilinos.org/docs/r11.12/packages/tpetra/doc/html/lesson01_no_mpi_8cpp-example.html
> [5]
>  [6]
> 
> http://trilinos.org/docs/r11.10/packages/tpetra/doc/html/classTpetra_1_1CrsMatrix.html#a139a91ad9430f012bfb3bcb8fe413ae1
> [6]
> 
> 
> 
> Links:
> ------
> [1]
> http://trilinos.org/docs/dev/packages/tpetra/doc/html/classTpetra_1_1CrsMatrix.html#a27c06cd71c7542742a731a6cc91180d2
> [2] http://www.feast-solver.org/
> [3] https://software.sandia.gov/mailman/listinfo/trilinos-users
> [4]
> http://trilinos.org/docs/r11.12/packages/tpetra/doc/html/classTpetra_1_1Map.html#aeb556362352e2ae62804601e36672a6c
> [5]
> http://trilinos.org/docs/r11.12/packages/tpetra/doc/html/lesson01_no_mpi_8cpp-example.html
> [6]
> http://trilinos.org/docs/r11.10/packages/tpetra/doc/html/classTpetra_1_1CrsMatrix.html#a139a91ad9430f012bfb3bcb8fe413ae1


More information about the Trilinos-Users mailing list