VecLoadIntoVector and double dispatch

Matthew Knepley knepley at
Fri Dec 4 08:57:17 CST 2009

On Fri, Dec 4, 2009 at 8:52 AM, Barry Smith <bsmith at> wrote:

> On Dec 4, 2009, at 2:58 AM, Jed Brown wrote:
>  On Thu, 3 Dec 2009 16:40:10 -0600, Barry Smith <bsmith at>
>> wrote:
>>>   The vector has a pointer to the DM so the VecView() for that
>>> derived vector class has access to the DM information.  The same
>>> viewer object can be used with a bunch of different sized Vecs since
>>> it gets the DM information out of each vector. With your model when
>>> the DM info is attached to the viewer you need to pass in a different
>>> viewer for different vectors (of different sizes). This is just plan
>>> nuts.
>> The DM info is not attached to the viewer, I have plain Vecs with a DM
>> composed so it can call DMView as in your snippet below.  The viewer
>> holds no state with regard to the DM.  But VecView_MyViewer_DM needs to
>> get called somehow,
>    Agreed.
>  it will in turn call DMView and use an internal
>> Viewer API so that the file has a connection from the Vec to the DM.
>>    I say this is nonsense, it ain't more closely coupled to the viewer.
>> VecGetArray is public, all vectors have to implement it in order to
>> function in PETSc.
>    This is not accurate. The SAMRAI vector class does not implement it.
> Yes, this means the SAMRAI vector class cannot use any PETSc built in matrix
> classes, but that is ok it provides its own. (note also that PETSc allows
> incomplete implementations of abstract classes where a new class
> implementaiton simply does not provide some methods and errors are generated
> if one tries to use them; I think this is a great design feature because it
> allows many-method base abstract classes, but I realize most of the
> computing world would be shocked and disgusted by this design feature).
>  The internal Viewer API to write the bytes of a Vec
>> with some connection to a DM is not public in the sense that users never
>> call it directly and no other viewers implement it.  You can substitute
>> any other Vec and this code will still write it correctly because your
>> Vec would be broken if it didn't implement VecGetArray [1].
>> DM currently monkey-patches Vec to set ops->view and
>> ops->loadintovector, but it doesn't rely on the implementation of Vec
>> (the definition of Vec_Seq/MPI is never used in src/dm) and in principle
>> DACreateGlobalVector could call VecSetFromOptions so that it's actually
>> operating with any Vec type (e.g. the Vec could implement it's
>> operations on a GPU, only bringing the values to the CPU for
>> VecGetArray).  Many viewers (formats) have no special representation of
>> vectors that come from a DM [2] which is the reason for having
>> ops->viewnative and ops->loadintovectornative.
>>  But somehow, calling VecView needs to be able to call the viewer
>>>> with enough context for the viewer to wire up the metadata.  Note
>>>> that the coordinates of the nodes in the DM are part of the vector
>>>> (it's an ALE formulation) so the metadata is nontrivial.
>>>   Each vector has a pointer to its DM (DMComposite etc) so its
>>> viewer can dump/read that info when needed. In fact, when writing a
>>> vector it could do something like
>>> myvecview(Vec v,Viewer viewer)
>>> {
>>>    DM dm = getreferencetoDMfrom(v);
>>>    DMView(dm,viewer);
>>>    ... now dump my vector with writes or whatever
>> My code looks exactly like this, my point is that the elided part
>> involves only the public Vec API, but a private Viewer API.  I can put
>> this chunk of code anywhere, but it needs changes when I make private
>> changes to the viewer, it doesn't need anything when I make private
>> changes to the DM or Vec (it doesn't even need to know the types of
>> those objects).  It is currently sitting with my DM, mostly because I
>> wrote the DMView first.
>>  You have to admit that it would be crazy to have the DM viewer method
>>> stored inside a viewer object which seems to be what your model is.
>> I wish I could get the dispatch with a symmetric relationship, but in
>> lieu of that, I don't especially care where the function pointer is,
>> what matters is that there is a way to get DMView_MyViewer_YourDM
>> called [3].  As long as dispatch is happening via PetscTypeCompare against
>> all the known types, we are restricted in that we either have no way to
>> use a third-party Viewer, or no way to view a third-party DM.  The
>> current choice is the former; I don't think it's inherently superior,
>> it's just a choice of which component can be extended without modifying
>> every object it touches.
>> If I was distributing a binary-only viewer plugin, I would have no
>> choice but to put all these methods in with the viewer (but I would
>> currently have no way to get them called).  Note that a viewer can be
>> very functional without ever using a private API or looking underneath
>> obj->data.  With all other PETSc objects (DM not so much because nobody
>> works with the generic DM interface), I can distribute a binary-only
>> implementation which is first class after dynamic loading.
>> My DMView is packaged up with my DM because it needs to know many
>> internal details of the DM, but the higher level logic about how to
>> write a time-dependent file with association between Vecs and DMs and
>> possibly other objects is currently with my viewer.  Now I need to write
>> DMView_MyViewer_DMComposite, but I won't be able to get it called
>> without modifying DMComposite (not okay since my viewer is not part of
>> PETSc).
>>    You know my answer. The Viewers are abstractions of file IO, ;
>>> except instead of having only methods like read()/write() bytes it has
>>> somewhat "higher level" operations. Note: these methods can include
>>> graphics operations for example that are still essentially like
>>> writes(). YOU are making the viewers into these strange "all knowing"
>>> beasts which I think is wrong.
>> Current PETSc viewers put minimal metadata and associations into a file,
>> the application is mostly responsible for "knowing" what was in the file
>> and wiring everything up correctly.  I'm trying to get more of this in
>> the file, partly to go with a VisIt plugin that is not tied to a
>> specific model.  I don't want to maintain a bunch of manual calls to
>> HDF5 to walk this file in the VisIt plugin especially when this
>> functionality is useful in other contexts [4].  I also hate C++.  ;-)
>> So I've been putting stuff like walking the header and attributes to
>> find out what's there, selecting a time step, etc., into my viewer.  The
>> viewer is able to offer up information like which fields are available
>> on which meshes, the plugin is a wrapper that uses this, DMLoadIntoDM
>> (but not this elegant), VecLoadIntoVector, etc., and stuffs it into
>> avt/vtk data structures.  Note that this plugin makes no direct calls to
>> HDF5 and uses only the public API of the Viewer, DM, and Vec.
>   Ok so you are deriving a new Viewer that implements much higher level
> operations that are particular to PETSc (since it knows about Vecs and DMs).
> This is perhaps reasonable. But the vecload/view methods still belong in the
> Vec object :-)
>> Jed
>> [1] In my view, a Vec is just an element of some finite dimensional linear
>> space.  If you want to use a special format, or stuff extra info in
>> there to make function evaluation easier or some operation faster,
>> that's fine, but it's still just an element of this linear space and it
>> better satisfy the appropriate axioms.  The linear space that this
>> vector lives in is defined by the DM and we need this to do anything
>> useful (like visualization) with the vector.
>   To the algebraic solvers I agree a " Vec is just an element of some
> finite dimensional ....." (since that is all the solvers need to know about
> them).  But to a specific application a Vec can be thought of as the linear
> algebra beast PLUS all kinds of grid information etc (like the DM).  So I
> think of a (application) Vec in the application as basically derived off of
> two classes:  the linear algebra beast PLUS a DM (it is implemented as a
> isA(linear algebra beast) and hasA(DM beast)).  This is why I feel so
> strongly that the VecView class that also reads/writes stuff related
> specifically to the Vec belongs with the Vec, not with the Viewer.

I would like to point out here that the notion of a Vec as a member of a
vector space is very often in adequate. What we really want is a section of
a fiber bundle. This
is much more appropriate for most computational settings. That is why we
want grid information to tag along with the Vec.


>  Regarding where code needs to be changed/where : The reason for the
> "native/super" support for the VecLoadIntoVector() (which I admit is
> terribly crude and ugly) is exactly so that one can implement a new
> loader/writer directly for a specific Vec/DM combination IN THE APPLICATION
> CODE or another library. One does not have to add it to src/vec/vec/impls
> Do you have suggestions on how one could improve the "native/super" thing
> to satisfy you?
> As I said before a Viewer in PETSc is just an abstraction/generalization of
> a file pointer with read/write methods on the file pointer. So it disturbs
> me that a method in a Viewer would even know what a Vec. It leads to these
> loops of knowledge that make code much more complicated. I would do anything
> to prevent these loops of knowledge I consider them so troublesome.
>   Barry
>> [2] Or the Vec doesn't use information about the DM when viewing on that
>> format if you like, I think this relationship is symmetric and
>> conceptually owned by neither.
>> [3] I think the location of the source for DMView_MyViewer_YourDM is
>> mostly a packaging and distribution matter rather than a design matter.
>> I don't see a problem with it being distributed independently of either
>> MyViewer or YourDM (except for the technical issue of how it announces
>> it's existence).
>> [4] It's all slightly complicated by the fact that for visualization,
>> it's important to be able to read chunks of the domain independently,
>> whereas the analysis code would prefer to only define proper collective
>> operations.

What most experimenters take for granted before they begin their experiments
is infinitely more interesting than any results to which their experiments
-- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <>

More information about the petsc-dev mailing list