← Back to team overview

dolfin team mailing list archive

Re: pointmap

 



Anders Logg wrote:
On Sun, Mar 25, 2007 at 04:46:49PM +0200, Garth N. Wells wrote:

Anders Logg wrote:
On Sun, Mar 25, 2007 at 10:44:58AM +0200, Garth N. Wells wrote:
Anders Logg wrote:
On Sat, Mar 24, 2007 at 06:52:04PM +0100, Garth N. Wells wrote:
In the old FFC output format, the function

   void pointmap(Point points[], unsigned int components[],
                 const AffineMap& map) const

returns in components[] a degree of freedom identifier (e.g. 0 for u, 1 for v, 2 for p, etc) for each entry in the element tensor. How can we get this information with the new UFC format? (or how can we avoid requiring it?)

Garth
Do you mean for evaluating dofs on user-defined functions (to get the
expansion coefficients in the nodal basis to put in the array w) or
for setting boundary conditions?

I think that in both cases it should be enough to evaluate the dofs on
the ufc::function, but we might have missed something. The function
evaluate_dof() takes a function f that may or may not be vector-valued
and computes the scalar value of dof i. So for a 2D vector-valued
Lagrange element of degree 1, dof 0 will be f_0(v0), dof 1 will be
f_0(v1), dof 2 will be f_0(v2), dof 3 will be f_1(v0) etc. The
function evaluate() in ufc::function needs to compute all values of
the possibly tensor-valued function.

What I don't see is how to make the link between a user-defined function in terms of x and j (j being u0, u1 and p for Stokes), and ufc::function which is in terms of x and i (i = 0 --> space_dimension).

Garth
The ufc::function needs to evaluate all values at once, not one at a
time. So for 2D Stokes, the array values would be of length 3.

We could keep the current interface eval(p, i) for Function and then
call this for each i to fill in the array values, but perhaps we
should change the interface of the class Function in DOLFIN so that
for a vector-valued function one needs to set all the values at once.
An advantage of this would be that one would not need the

   if ( i == 0 )
     return 3.0;
   else if ( i == 1 )
     return 5.0;
   else
     return 6.0;

Instead, one would write

   values[0] = 3.0;
   values[1] = 5.0;
   values[2] = 6.0;

The array values is a contiguous array that may represent something
tensor-valued as well. If the function takes values in R^{3x3} then
the array has length 9.

For a vector and scalar function this would be be nicer, but we have to make sure that it works nicely with mixed elements. When applying boundary conditions, not all vector components will always be supplied and something like an boolean array will be required to denote Dirichlet bc's.

This is more difficult than I thought... The problem seems to be that
if we use evaluate_dof from the UFC interface to get the value of dof
number i located on the boundary when evaluated on a given function f
(the Dirichlet boundary condition), we want to have the option of not
saying anything about the value for certain *components* of the
function (for which we want Neumann or do nothing boundary
conditions). So we would need a mapping from component to dof in order
to do this, but we don't have this and I don't think we can have it in
general. (If the element is not a simple tensor product of scalar
elements like vector Lagrange.)


The vast majority of analyses use Lagrange basis functions, so I think that we should support it, and for more complex basis functions the most flexible approach would be weak enforcement of Dirichlet boundary conditions.

Perhaps we could solve the problem by requiring that a boundary
condition must always be given by a Function? And that the boundary
condition applies to all components that the Function represents. (No
option of not setting the boundary condition.)

Sounds like something as simple as a slip boundary condition would be become complex (or impossible) to impose.

In the case when one
has a mixed system, then one needs to supply one Function for each sub
system that one wants to set Dirichlet conditions for, and don't set
conditions for sub systems that shouldn't have a Dirichlet
condition. So one could do something like

  set_bc(A, x, u); // u is a Function
  set_bc(A, x, p); // p is a Function

This still does not solve the problem when one wants to set a
Dirichlet condition for u in one part of the domain but not in
another, but then we should support setting boundary conditions for
sub domains which are specified by a MeshFunction (over facets).


I would like to see the current functionality as well as the possibility to set boundary conditions on subdomains (very useful when the geometry is not exact).

Garth


On a related note, I have mentioned earlier that Kalkulo (kalkulo.no)
is developing a mesh tool for DOLFIN (paid for by Simula). They are
still working out a few bugs, but it should soon be ready (free and
available from fenics.org as part of FEniCS). The mesh tool lets you
mark (graphically) subsets of the boundary and writes the results as a
DOLFIN XML MeshFunction file.

/Anders
_______________________________________________
DOLFIN-dev mailing list
DOLFIN-dev@xxxxxxxxxx
http://www.fenics.org/mailman/listinfo/dolfin-dev




Follow ups

References