Subject: [FFC-dev] Mesh

From: <[email-masked]>
Date: Wed, 31 Aug 2005 05:53:08 -1100

I'm afraid that this discussion of what to name various pieces
of information about the mesh will lead to code growth by
accretion-on-demand rather than design.  It will also lead us
into maintaining separate code branches in 2d and 3d. 
Consider instead the following:

First, every mesh entity, whether a vertex, an edge, etc, has
at least two intrinsic geometric features:
1.) its intrinsic topological dimension (a vertex is dimension
0, and so on)
2.) its measure in that dimension (length of an edge, area of
a triangle, etc)

This suggests that "vertex","edge", etc are really instances
(or subclasses if you must) of an abstract "mesh entity" class
that have as member information the topological dimension and
coordinates.

Moreover, every mesh entity has certain topological/connection
information associated with it.  We might need to know what
faces are attached to which tets, which tets share a given
face (think discontinuous Galerkin or high degree Lagrange),
etc.  This connection information can be stored in a single
dictionary for each mesh entity -- Matt Knepley calls it an
"incidence relation".  This dictionary maps topological
dimension to a list of entities of that dimension that are
connected to the given entity.  So, for a vertex, you put "1"
into the indicence dictionary and get a list of all the edges
for which that vertex is an endpoint.  You put "2" into the
dictionary and you get out a list of triangles that form a
patch around the vertex.  For a triangle, you can put in "2"
and get what triangles share an edge, and so on.

This gives a minimal interface that you don't have to maintain
separately for different dimensions of meshes, and you don't
have to haggle over names.  All you have to do is get the
information from the mesh generator and into these data
structures.  Matt's package has this all done (plus
parallelism).  You should keep an eye on when this will be
ready to plug in.

Rob

---- Original message ----
>Date: Wed, 31 Aug 2005 17:57:50 +0200 (MEST)
>From: [email-masked]  
>Subject: Re: [FFC-dev] Tensors  
>To: "Discussion of FFC development" <[email-masked]>
>
>> I'm planning to add a list of predefined operators and
operands to FFC
>> to simplify and beautify the notation anyway, so this would
fit well
>> into that category.
>>
>> Operators that will be added include grad(), rot() and div().
>>
>> Operands that will be added include h (mesh size) and d
(space dimension).
>>
>> Any other suggestions?
>
>I find it useful to be able use some more geometry
information (in
>particular when we start integrating on surfaces):
>
>vol - volume of cell
>area(i) - area for local face i
>hface(i) - h for local face i
>...
>same for local edges?
>
>It may also be useful to be able to identify a local face from a
>collection of vertices; that is, while visiting a cell in the
assembly one
>might want to get information about the local faces.
>
>On the other hand, we may avoid a lot of this since we intend
to assemble
>separately over faces (and edges), and then vol and h may be
reduced to
>the corresponding functions for the face when looping over
the faces.
>
>/Johan
>
>
>
>
>
>
>
>_______________________________________________
>FFC-dev mailing list
>[email-masked]
>http://www.fenics.org/cgi-bin/mailman/listinfo/ffc-dev

---

From: Anders Logg <[email-masked]>
Date: Wed, 31 Aug 2005 15:16:34 -0500

Yes, the intention is that we should look closely at Matt's new mesh
component and if possible completely throw out all our mesh stuff and
wrap the new mesh component in DOLFIN.

When you're ready Matt, pick a name (Sieve?) and put it on fenics.org.

My only worries is that (i) the new mesh component will be not be
limited enough in its scope and force us to throw out more than just
the mesh component and (ii) it will be impossible to install... ;-)

/Anders

-- 
Anders Logg
Research Assistant Professor
Toyota Technological Institute at Chicago
http://www.tti-c.org/logg/

---

From: [email-masked]
Date: Thu, 1 Sep 2005 08:22:43 +0200 (MEST)

Rob,

All this is great, and it should be interesting to see Matt's mesh
package. I'm sure it will be very useful for FEniCS. My understanding of
our current discussion is that it was related to what mesh information
that should be avaliable through FFC when defining a form. Not a
discussion of the possible names for this functionality, I certainly do
not suggest to include functions like "h_face_x" or similar. I think the
question was what information about the mesh needed to be avaliable
through FFC.

/Johan

---

From: Matthew Knepley <[email-masked]>
Date: Thu, 01 Sep 2005 09:51:37 -0500

> When you're ready Matt, pick a name (Sieve?) and put it on fenics.org.
  
  There will be two names: Sieve for a thing which encodes topology, and
SievedArray for a thing which stores data in a hierarchical format consistent
with a certain topology. The idea here being that a SievedArray can both

  1) Present data in a contiguous parallel format like a PETSc Vec

and

  2) Restrict that data to coverings of the original topology (like a
     single element), and then reassemble

> the mesh component and (ii) it will be impossible to install... ;-)

  I will try and address these:

  1) I am fairly confident that we are alright here. This is ONLY intended
     to make the program understand topology and data organized around it.
     It should slot in just below normal linear algebra for storage, and
     should provide the necessary information to organize a) parallel
     assembly and b) finite element assembly.

  2) Since it is implemented entirely in C++ with STL as a first pass, our
     idea was for you just to use the C++, which is a bunch of self-contained
     files. The build should be easy. I have a bunch of interop stuff on top
     to make it do Python, but you don't need this.

     Matt

-- 
"Failure has a thousand explanations. Success doesn't need one" -- Sir Alec Guiness

---

From: Anders Logg <[email-masked]>
Date: Thu, 1 Sep 2005 10:36:05 -0500

This really sounds good. I look forward to trying it out.

The array part seems to correspond to the Function concept we have in
DOLFIN. A Function in DOLFIN is defined by three things (private
member of the Function class):

    Vector* x; // wrapper for PETSc Vec
    FiniteElement* element;
    Mesh* mesh;

The vector holds the degrees of freedom, and the finite element
(generated from a FIAT finite element by FFC) defines how these are
distributed on the mesh.

A Function can be evaluated (but the implementation for evaluation at
arbitrary points is missing), passed as a argument to a variational
form and be written to a File for visualization.

/Anders

---

From: Matthew Knepley <[email-masked]>
Date: Thu, 01 Sep 2005 10:53:37 -0500

> form and be written to a File for visualization.

  Yes, the generalization comes from the broader definition of restriction
and assembly operations. In the language of the paper, a Function can define
only cocycles, where a SievedArray can handle cochains.

     Matt

---

From: Robert C. Kirby <[email-masked]>
Date: Fri, 2 Sep 2005 08:35:29 -0500

Johan, thanks for the clarification, but I stand by my point with FFC.  
I fear we are "accreting" onto it.  My complaints about C++ are 
organized on the same lines -- it gives the feel of having lots of 
stuff added to it as Stroustrup decided to do it over many years.  You 
get weird syntax, clashing of concepts, and it creates so much 
cognitive dissonance for me.

The point of the mesh package is that we should strive to make our mesh 
as well as what goes on in the FFC syntax well-designed, general, and 
elegant, rather than just saying, "I need to solve this problem, so 
I'll slap on some incremental changes for now".  While things may be 
implemented in that way (we have finite time), the design process 
should be more careful.  Some weak forms, such as stabilizing terms, 
refer to the mesh (such as "h"), and we should not train FFC to worry 
what dimension things are in its syntax any more than is absolutely 
necessary.

Robert C. Kirby
Assistant Professor
Department of Computer Science
The University of Chicago
http://people.cs.uchicago.edu/~kirby

---

From: Anders Logg <[email-masked]>
Date: Fri, 2 Sep 2005 09:06:16 -0500

Good point, but all things that I'm planning to add onto FFC will just
be shortcuts to things that one can already do today.

One example is the mesh size h. One can do stabilization already today:
Garth has implemented stabilized convection-diffusion in DOLFIN with
FFC, and Johan has stabilized Navier-Stokes. To do this, they just
declare a piecewise constant Function h in FFC (which is part of the
basic functionality, to have arbitrary user-defined functions in the
forms). The problem is that they then need to create a Function h in
DOLFIN and manually implement that the value of h should be the mesh
size on each cell. To simplify, I just need to add a predefined
Function in FFC, something like

    h = Function("mesh size")

This won't add to the complexity of FFC, since it's just an add-on
that works on top of the existing machinery.

Same thing with adding new operators grad, div, rot. They can also be
defined in terms of the existing language.

/Anders

---

From: Robert C. Kirby <[email-masked]>
Date: Fri, 2 Sep 2005 09:12:25 -0500

Correct I don't think stabilization presents a particular problem to 
FFC (other than generating a nonlinear system if you make it dependent 
on the residual) -- but I am making precautions as much as particular 
suggestions here.

---

From: [email-masked]
Date: Fri, 2 Sep 2005 16:55:20 +0200 (MEST)

The same thing is true for operations like 'sqrt' and division. Today we
declare a Function in DOLFIN and manually implement the inverse of the
magnitude of a certain Function w, which is sort of messy. It would be
much cleaner (and be more memory efficient) to be able to write
'1.0/sqrt(sqr(w1) + sqr(w2))' in FFC directly.

I agree that we should be very careful about the design. But that does not
have to be in conflict with having some nice functionality on top of that.

/Johan