dolfin team mailing list archive

[Jed Brown <jed@xxxxxxxx>]

On Mon 2008-08-25 08:18, Anders Logg wrote:
> On Sun, Aug 24, 2008 at 08:16:59PM -0400, Shawn Walker wrote:
> >
> > On Sun, 24 Aug 2008, Anders Logg wrote:
> >
> >> On Sun, Aug 24, 2008 at 03:55:49PM -0400, Shawn Walker wrote:
> >>>
> >>> On Sun, 24 Aug 2008, Anders Logg wrote:
> >>>
> >>>> The main problem here seems to be how to even represent a higher-order
> >>>> mesh. Are there any existing conventions for this?
> >>>
> >>> Not that I know of, other than the 'obvious' way.
> >>>
> >>>> Doing it locally is fine. For each cell, we have a convention for how
> >>>> we number vertices, edges and faces (it's in the UFC manual) but we
> >>>> don't have a convention for how to number mesh entities globally.
> >>>>
> >>>> We could use the dof map (for Lagrange elements) that FFC computes
> >>>> based on the numbering of mesh entities, but that computation is based
> >>>> on a global numbering of edges and faces, and that numbering depends
> >>>> on the algorithm for computing these entities (in TopologyComputation).
> >>>>
> >>>> So, the only solution I see is to store the coordinates cell-wise,
> >>>> which would mean duplicating the data, both in the XML file and for
> >>>> storing the coordinates in MeshGeometry.
> >>>
> >>> I think this would be ok.  But it must still be possible to update the
> >>> mesh geometry (i.e. update the coordinates of all the control points for
> >>> the curved triangle geometries) in, say, a time-stepping loop.  For
> >>> example, you would need this for ALE moving mesh applications.  This
> >>> means you must update the individual cell geometries in a 'compatible'
> >>> way.  I haven't actually thought about this.  Would this be do-able?  If
> >>> so, I have no problem with keeping the higher order mesh geometry
> >>> completely local.  So the question is: "Is it possible to access and
> >>> change the local higher order mesh geometry such that the higher order
> >>> mesh coordinates stay consistent?"
> >>>
> >>> - Shawn
> >>
> >> Yes, we just need to have a global list of coordinates, and then
> >> mappings from all cells to the coordinate indices of those cells.
> >>
> >> Then one can move a coordinate once and it will affect all cells
> >> sharing that coordinate.
> >>
> >> My suggestion would be to have one single list of coordinates like we
> >> do now in MeshGeometry. The first N coordinates would be the vertex
> >> coordinates (like now) and then the remaining coordinates (if any)
> >> would be stored per cell, so that coordinate N is the coordinate for
> >> the first edge midpoint of the first cell, coordinate N + 1 is the
> >> coordinate for the second edge midpoint of the first cell etc.
> >
> > ok, just to be clear, you mean:
> >
> > Have a list of coordinates of length N+Q, where 0 through N-1 is for the 
> > coordinates used now, and the remaining N to N+Q-1 coordinates contains 
> > coords that lie on the faces/edges of the mesh.  And N through N+Q-1 is 
> > ordered by cell index grouping.
> >
> > That seems fine.  How would the .xml file look?  Perhaps this:
> >
> > <?xml version="1.0" encoding="UTF-8"?>
> >
> > <dolfin xmlns:dolfin="http://www.fenics.org/dolfin/";>
> >   <mesh celltype="triangle" dim="2" map_type="P2">
> >
> >     <vertices size="4">
> >       <vertex index="0" x="0.0"  y="0.0"/>
> >       <vertex index="1" x="1.0"  y="0.0"/>
> >       <vertex index="2" x="0.0"  y="1.0"/>
> >       <vertex index="3" x="-1.0" y="0.0"/>
> >     </vertices>
> >
> >     <cells size="2">
> >       <triangle index="0" v0="0" v1="1" v2="2"/>
> >       <triangle index="1" v0="3" v1="0" v2="2"/>
> >     </cells>
> >
> >     <higher_order_vertices size="5">
> >       <vertex index="4" x="0.5"  y="0.05"/>
> >       <vertex index="5" x="0.6"  y="0.6"/>
> >       <vertex index="6" x="0.0"  y="0.5"/>
> >       <vertex index="7" x="-0.5" y="0.5"/>
> >       <vertex index="8" x="-0.5" y="0.0"/>
> >     </higher_order_vertices>
> >
> >     <higher_order_cells size="2">
> >       <triangle index="0" affine="false" v0="4" v1="5" v2="6"/>
> >       <triangle index="1" affine="true"  v0="8" v1="6" v2="7"/>
> >     </higher_order_cells>
> >
> >   </mesh>
> > </dolfin>
> >
> > Or should it be grouped by triangle face?  This also should be  
> > extensible to say cubic edges, or even higher.
> 
> Yes, this looks good. But maybe the local indices for the higher-order
> cells should not start at zero:
> 
>      <higher_order_cells size="2">
>        <triangle index="0" affine="false" v3="4" v4="5" v5="6"/>
>        <triangle index="1" affine="true"  v3="8" v4="6" v5="7"/>
>      </higher_order_cells>

For what it's worth, mesh formats like exodusII and VTK accept quadratic
elements.  They view a quadratic element as a different topology, for
instance there is VTK_HEXAHEDRON=12 and VTK_QUADRATIC_HEXAHEDRON=25.  Of
course each of these has a canonical numbering for the nodes.  All
vertices are kept in one array, then connectivity for each topology
(either a separate array for each topology or interlaced with a topology
ID for every element).

I'm not aware of any standard for cubic and higher order meshes (or any
visualization software that could use it).  In that case, I think the
most reasonable thing is to store full element connectivity without
reference to vertices, that is (Region->Face), (Face->Edge),
(Edge->Vertex) and then store coordinates as a function over this space
(however you store such things, such as by indexing interior degrees of
freedom for every element/facet).  Note that there are many more
coordinates than vertices.

Providing a canonical numbering for nodes on say, a cubic hex, requires
98 vertices.  Also, the mesh function approach plays nicely with mixed
order discretization, modal bases, and nonconforming meshes.

So it might be reasonable to just implement the quadratic topologies now
and wait for cubic and higher order since it is a major undertaking to
do in a flexible manner.

Jed

Attachment: pgptTDIn45ShA.pgp
Description: PGP signature