The biggest reason we haven't added "multi-vectors" is more design
decisions that have to be made then actual code implementation. There
are two important design decisions:
I) Introduce a new class or overload Vec?
1) Introduce a Vecs class, essentially copying the interface from
Vec and implement the methods. A Vecs is essentially a collection of Vec.
Then provide methods such as
MatMults(), MatSolves(), PCApplys(), PCSetUps(), KSPSolves() etc
2) Overload the current Vec class so that a Vec is either a "regular"
Vec or a collection of vectors. Modify the Vec interface and methods
as needed, then modify the MatMult() etc as needed
3) Do 1) and have it all working well and then use this to quickly
do 2) and throw away the seperate concept of Vecs (without losing
the fastest possible performance for a "new" Vec that has a single
vector).
The advantage of 1) is we never break the current interface. 2) has the
advantage of not introducing a new class; and since software
comprehension difficulty grows exponentially with the number of classes
this is a good thing. The advantage of 3) is we get an incremental approach
to 2) (incremental is often good). Disadvantage of 3) is that users have
to live with first the Vecs interface (but on the other hand most people
won't use Vecs anyways) and then go back to Vec.
I think 3) is the correct approach.
II) Default storage of "multivectors"? (note: the interfaces for Vecs will, of
course, be data storage format neutral, but the code we actually write
cannot be)
1) Have the Vecs be essentially an array of Vec.
2) Have a Vecs essentially be a double array where the "values" of the
"vectors" are interlaced.
The advantage of 1) is that it is easier to implement and doing things
like deflation (removing one of the vectors in a Vec) doesn't require
any copies of the vector data. The advantage of 2) is better performance,
and since we are using block for better performance we should take advantage
of the best possible performance.
I think implementing both 1) and 2) is appropriate, starting with 1) since
getting it completed would be faster.
Barry
Note: the Vecs currently in petscvec.h is just a toy and can be ignored.
On Fri, 7 Apr 2006, Victor Eijkhout wrote:
> The question has come up before. If all goes according to plan, then Bill
> Gropp will hire a student for the summer to implement it. On a project that I
> lead, and that has use for multivectors for a different purpose, there is
> actually budget. Feel free to recommend someone to Bill if you know of a good
> candidate. Does your former advisor have students who are familiar with the
> matter?
>
> Victor.
>
> On Apr 7, 2006, at 1:05 PM, Richard Tran Mills wrote:
>
>> Hi Guys,
>>
>> My former Ph.D. advisor, Andreas Stathopoulos (andreas(a)cs.wm.edu) has
>> released an eigensolver package, PRIMME, that he is interested in
>> interfacing as an external package with PETSc. I told him a bit about how
>> he could do this, but what I couldn't really answer was his questions about
>> "multivector" support in PETSc. PRIMME implements block
>> Davidson/Jacobi-Davidson type methods, and thus needs to do a lot of
>> operations like sparse matrix-multivector multiplies. I haven't seen
>> anything like this in PETSc (and since PETSc doesn't have block Krylov
>> solvers, I suppose there is currently no need). Has anyone done any work
>> on stuff like this, or does anyone have any plans to? Any thoughts on
>> this?
>>
>> Cheers,
>> Richard
>>
>
>
