[petsc-dev] Fwd: DMplex reader / viewers

[Garth N. Wells]

On 2014-01-22 11:32, Gorman, Gerard J wrote:
> Hi Garth
> 
> I was wondering what your view was on this whole discussion since you
> and Richard integrated XDMF with Dolfin in a dCSE last year
> (http://www.hector.ac.uk/cse/distributedcse/reports/UniDOLFIN/UniDOLFIN/index.html
> [2]). With the benefit of hindsight are you happy with the approach
> you took? Was the XDMF metadata expressive enough?

XDMF is great for post-processing (it's supported by ParaView and 
VisIt), and it's good for reading in vanilla meshes (i.e., no markers). 
It's much, much better than the native VTK formats for visualisation - 
it's two files only and for time dependent problems the mesh only needs 
to be written once. There are some issues with how it orders points - 
the index of a vertex is implicitly taken as its position in the file. 
This isn't ideal if a user wants to inspect a vertex or cell with a 
particular index in ParaView. We make sure we order correctly in the 
file and we remove duplicates on process boundaries to avoid 
visualisation artefacts. With the HDF5 storage, we get orders of 
magnitude speeds-ups over VTK/XML for IO.

It might be possible to add markers, etc to XDMF in a consistent way. We 
haven't looked at that.

> It seems like you
> must have ran into some problems as IO for vectors was not supported
> for example.
> 

For IO of other objects, we bake our own HDF5 files. It would prefer to 
work against a specification, but that area is a big mess . . .

Garth

> Cheers
> Gerard
> 
> Begin forwarded message:
> 
>> FROM: "Gorman, Gerard J" <g.gorman at imperial.ac.uk>
>> 
>> SUBJECT: RE: [PETSC-DEV] DMPLEX READER / VIEWERS
>> 
>> DATE: 21 January 2014 23:58:11 GMT
>> 
>> TO: "Tim Tautges (ANL)" <tautges at mcs.anl.gov>
>> 
>> CC: petsc-dev <petsc-dev at mcs.anl.gov>
>> 
>> On 21 Jan 2014, at 21:07, Tim Tautges (ANL) <tautges at mcs.anl.gov>
>> wrote:
>> 
>> On 01/21/2014 12:00 PM, petsc-dev-request at mcs.anl.gov wrote:
>> 
>> Date: Tue, 21 Jan 2014 17:56:35 +0000
>> From: Blaise A Bourdin <bourdin at lsu.edu>
>> To: "Gorman, Gerard J" <g.gorman at imperial.ac.uk>
>> Cc: petsc-dev <petsc-dev at mcs.anl.gov>
>> Subject: Re: [petsc-dev] DMplex reader / viewers
>> Message-ID: <ECA3EBD3-7618-490A-B167-B562428ECC7B at lsu.edu>
>> Content-Type: text/plain; charset="Windows-1252"
>> 
>> On Jan 21, 2014, at 10:23 AM, Gorman, Gerard J
>> <g.gorman at imperial.ac.uk> wrote:
>> 
>> Hi Blaise
>> 
>> This is pretty much music to my ears. Interoperability is a major
>> headache, most file formats are rubbish (e.g. cannot support for
>> boundaries, limited meta data etc etc) most lack reasonable parallel
>> IO solutions (i.e. something other than 1 file per process)?and so
>> on. At the moment a few groups here at IC (this mostly means Michael
>> Lange ;) are integrating DMPlex with a range of FEM codes to try to
>> rationalise this. There are certainly some features missing (we?ve
>> bumped into most of that you?ve listed below) but all indications
>> are that DMPlex has the right design to be expanded to cover these
>> use cases.
>> 
>> On 21 Jan 2014, at 15:30, Blaise A Bourdin <bourdin at lsu.edu> wrote:
>> 
>> Hi,
>> 
>> It looks like DMplex is steadily gaining maturity but I/O is
>> lagging behind. As far as I understand, right now, PETSc can _read_
>> a mesh in exodus format, and write binary VTS format, but many
>> issues remain, IMHO:
>> - The exodus reader relies on a hard-coded nodeset named ?marker?.
>> Generating such a nodeset is not trivial
>> (at least not for complex meshes generated with Cubit / Trelis).
>> - Reading from or writing to exodus files is not supported.
>> 
>> So we *really* want this as well for the purpose of dumping results
>> and checkpointing. Trying to shoehorn high order and DG data into
>> VTK is complete crap. I love VTK because its got loads of
>> functionality and it is easy to throw stuff together, but we are
>> getting ourselves into the position that we are just layering hack
>> after hack to deal with the fact that we cannot write all required
>> data into a vtk file. These days VTK/paraview has its own ExodusII
>> reader so we have a route to nice seamless integration.
>> 
>> - The VTS viewer only allows to read and write _all_ fields in a
>> DM. This may be overkill if one only
>> wants to read boundary values, for instance.
>> 
>> or only writing out prognostic fields for example.
>> 
>> - The VTS viewer loses all informations on exodus nodesets and cell
>> sets. These may have some significance
>> and may be required to exploit the output of a computations.
>> 
>> Right - this includes boundary labels, and it is just a fluff to
>> have to write this out into VTK. You would have to write a separate
>> vtu or something resulting in more messiness (and I already have
>> enough problems on LUSTER from having too many files).
>> 
>> - VTS seems to have a concept of ?blocks?. My understanding is that
>> the parallel VTS viewer uses blocks to
>> save subdomains, and that continuity of piecewise linear fields
>> across subdomain boundaries is lost.
>> It is not entirely clear to me if with this layout, it would be
>> possible to reopen a file with a
>> different processor count.
>> 
>> I think you just do not want to go there? For a start your vtk file
>> would not be a valid checkpoint to consider restarting on a
>> different number of processes. And secondly, it would just be a mess
>> to program.
>> 
>> I can dedicate some resources to improving DMplex I/O. Perhaps we
>> can start a discussion by listing the desired features such readers
>> / writers should have. I will pitch in by listing what matters to
>> me:
>> 
>> Keep talking?we have also an FTE working on this currently but this
>> is a long wish list and a lot of effort is required if this is to be
>> done within a reasonable time frame. It would be great if more
>> people were working on this.
>  I have a postdoc here who will devote some of his time to this task.
> 
>>> - A well documented and adopted file format that most
>>> post-processors / visualization tools can use
>> 
>> ExodusII appears to be the current favoured option.
>  Sadly yes? But SILO may be a close second and has a more modern
> interface.
> 
>>> - Ability to read / write individual fields
>>> - Preserve _all_ information from the exodus file (node / side /
>>> cell sets), do not lose continuity of fields
>>> across subdomain boundaries.
>>> - Ability to reopen file on a different cpu count
>> 
>> So this is where we need to have parallel IO support. Quoting from
>> petcs?s exodus.py
>> ??"
>> # ExodusII does not call HDF5 directly, but it does call
>> nc_def_var_deflate(), which is only
>> # part of libnetcdf when built using --enable-netcdf-4. Currently
>> --download-netcdf (netcdf.py)
>> # sets --enable-netcdf-4 only when HDF5 is enabled.
>> ??"
>> So, there may be some rebuilding required to ensure that all the
>> dependencies are built properly but it?s clearly there.
> 
>  I am not sure if Exodus has a good solution here. As far as I
> understand, exodus is inherently sequential, even when implemented
> with HDF5 instead of netcdf. I would also worry about third party
> support for exodus files using HDF5 as their storage format.
>  Exodus has an parallel extension called nemesis, but I can?t figure
> out how how their concept of ghost point and cells works. The
> documentation on this point is really unclear.
> 
>  I have to admit I was kind of hoping that ExodusII folks would have
> come on a bit more on the parallel IO front (I’m assuming those guys
> also run large simulations…). That said, I see this as a two stage
> process: first integrate with DMPlex as that should give the high
> level abstraction for read/write to file; secondly extend the family
> of readers/writers. At least this way there will be some agility and
> interoperability between different formats, and it will not be too
> disruptive to the application codes when a different formats adopted.
> 
>>> Considering the dismal state of parallel FE formats and libraries,
>>> it seems to me that one needs to chose between two options:
>>> 
>>> a. scatter back to a single I/O node and use sequential I/O using
>>> the ordering of the original (exodus) mesh. This allows reading
>>> and writing on an arbitrary number of processors, but has
>>> potential memory footprint and performance issues. How large a
>>> mesh can we reasonably expect to be able to handle this way?
> 
>  Personally I would stay far away from this option. Other than being a
> terrible serial bottleneck, it’s a major headache when you just want
> to run something just a little bit bigger than what happens to fit
> within a single node…
> 
>>> b. Do ?poor man? parallel I/O where each CPU does its own I/O, and
>>> possibly create interface matching files ? la nemesis or SILO.
>>> Maybe, we can save enough information on the parallel layout in
>>> order to easily write an un-partitionner as a post-processor.
> 
>  I am pretty sure that if we are writing everything in slabs to a HDF5
> container we do not have to worry too much about the parallel layout
> although some clear optimisations are possible. In the worst case it
> is a three stage process of where we perform a parallel read of the
> connectivity, scatter/gather for continuous numbering, parallel
> repartitioning and subsequent parallel read of remaining data.
> Importantly, it is at least scalable.
> 
>>> Unless one can come up with a better solution than a or b, I?d
>>> like to start by writing a very simple ASCII viewer demonstrating
>>> the communication involved, then modify it to use exodus, SILO or
>>> HDF5 format.
>>> 
>>> Blaise
> 
>>> 
>> 
>> Depending on the degree of direct interaction/automation in those
>> interactions between the mesh and Petsc, there are other options as
>> well. One that we're developing, based on the MOAB library, can
>> read/write (in serial) ExodusII, and also supports parallel
>> read/write using its own HDF5-based format. Parallel I/O robustness
>> has been iffy above ~16k procs and 32M-64M hex/tet elements, but for
>> smaller problems it should work. We're in the process of developing
>> direct support for going between a mesh defined with fields (called
>> tags in MOAB) and petsc vectors. MOAB has pretty solid support for
>> things like computing sharing and ghosting between procs and
>> exchanging/reducing field values on those entities. Viz is supported
>> either by compiling a VTK/Paraview plugin that pulls the mesh/fields
>> through MOAB or by translating to VTK (also supported directly from
>> MOAB); Visit also has a plugin you can enable. See
>> http://trac.mcs.anl.gov/projects/ITAPS/wiki/MOAB [1] for details of
>> MOAB; the petsc integration stuff is on a bitbucket branch of petsc
>> owned by Vijay Mahadevan.
> 
>  Another reason this could be of great interest is that MOAB supports
> (according to the docs) geometric topology which could be used when
> adapting the mesh on a curved surface for example - another item on my
> which list. Is it integrated to PETSc via the plex or does this
> essentially replace the functionality of the plex? Why does it break
> down for more than 16k procs? is it just a case that Lustre gets
> hammered? What magic sauce is used by high order FEM codes such as
> nek500 that can run on ~1m cores?
> 
>> libmesh also maintains its own DMlibmesh, but I'm not sure how
>> solid their support for large mesh / parallel I/O is (but they've
>> been working on it recently I know).
> 
>  Are there any other formats that we should be considering? It’s a
> few years since I tried playing about with CGNS - at the time its
> parallel IO was non-existent and I have not seen it being pushed
> since. XDMF looks interesting as it is essentially some xml metadata
> and a HDF5 bucket. Is anyone championing this?
> 
>  Cheers
>  Gerard
> 
> 
> Links:
> ------
> [1] http://trac.mcs.anl.gov/projects/ITAPS/wiki/MOAB
> [2]
> http://www.hector.ac.uk/cse/distributedcse/reports/UniDOLFIN/UniDOLFIN/index.html