dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

I agree and this is what I suggest. FFC generates a function
evaluate_dof (part of the UFC specification) which is exactly your
\int_\Gamma \eta_i if I understand it correctly. So, evaluate dof i on
g which is

    \int_\Gamma \eta_i(g) = evaluate_dof(i, g)

and then set u_i equal to this value.

/Anders


On Sun, Apr 08, 2007 at 09:34:21AM -0500, Matthew Knepley wrote:
> So I have been thinking about this myself, meaning how to apply Dirichlet
> conditions. I think the correct way comes from the FIAT representation, bu
> you guys can tell me what you think. I say something like
> 
>   u |_\Gamma = g
> 
>   \sum_j u_j \phi_j |_\Gamma = g
> 
> and we have a nodal basis, so we can use the dual basis eta_i to get
> 
>   u_i = \int_\Gamma \eta_i g
> 
> I believe this should work for any element. This does require quite a bit
> of code generation, but I did not think we were shy about that :)
> 
>   Matt
>
> On 4/8/07, Anders Logg <logg@xxxxxxxxx> wrote:
> > On Sun, Apr 08, 2007 at 11:27:47AM +0200, Garth N. Wells wrote:
> > >
> > >
> > > Anders Logg wrote:
> > > > On Sat, Apr 07, 2007 at 02:35:07PM +0200, Garth N. Wells wrote:
> > > >
> > > >> Any further thoughts on how to apply Dirichlet boundary conditions? I've
> > > >> been playing around with it, and with the new UFC interface I don't see
> > > >> how we can retain the simplicity of the old approach.
> > > >>
> > > >> Garth
> > > >
> > > > How about the following:
> > > >
> > > > - A boundary condition can be applied in one of two ways: either by
> > > > overloading eval() in the class BoundaryCondition, or by suppliying a
> > > > Function that defines the boundary condition.
> > > >
> > > > - In the first case, the eval() function looks the same as the eval()
> > > > method in Function, but instead of real we have BoundaryValue:
> > > >
> > > >     void eval(BoundaryValue* values, const real* x);
> > > >
> > > > When using this method, one sets the values of dofs direcly, so for a
> > > > BDM element, one would specify the value of the normal component, not
> > > > the function value.
> > > >
> > >
> > > Don't we need a function like pointmap to do this? tabulate_coordinates
> > > will give the position, but we are still missing something that
> > > translates the component of a vector function (0, 1 or 2) to the
> > > appropriate dof (like pointmap did).
> >
> > Ah, yes. I forgot about this already... Since it's not always possible
> > to create such a list (only for Lagrange), I would like to avoid it
> > (and preferrably also tabulate_coordinates) and just use evaluate_dof.
> >
> > If you look in the generated code for PoissonSystem.form (one of the
> > FFC demos), we have the following code in evaluate_dof:
> >
> >   virtual double evaluate_dof(unsigned int i,
> >                               const ufc::function& f,
> >                               const ufc::cell& c) const
> >   {
> >     ...
> >
> >     // Nodal coordinates on reference cell
> >     static double X[6][2] = {{0, 0}, {1, 0}, {0, 1}, {0, 0}, {1, 0}, {0, 1}};
> >
> >     // Components for each dof
> >     static unsigned int components[6] = {0, 0, 0, 1, 1, 1};
> >
> >     ...
> >
> >     // Evaluate function at coordinates
> >     f.evaluate(values, coordinates, c);
> >
> >     // Pick component for evaluation
> >     return values[components[i]];
> >   }
> >
> > FFC knows that the element is a tensor-product Lagrange element so it
> > generates the components array, evaluates the function and picks the
> > dof from the appropriate component. It would be good if we didn't need
> > to worry about the components in DOLFIN.
> >
> > So then my suggestion would be to skip the first option and require
> > that one always sets a boundary condition by a Function. This can be
> > very flexible. The Function can either be something user-defined (with
> > an overloaded eval) or a DiscreteFunction obtained as the solution
> > from a previous problem.
> >
> > The problem is then that we need to be able to set boundary conditions
> > only on a certain sub domain of the mesh, but this should be right way
> > to do it anyway. The current coordinate-based approach (if x == 0
> > then...) has only been a temporary solution until we figured out how
> > to label sub domains (by using a MeshFunction<uint>).
> >
> > Perhaps we could provide some functionality for defining sub domains
> > based on coordinates?
> >
> > /Anders
> >
> >
> > > > - The class BoundaryValue remains the same as today. So if one does
> > > >
> > > >     values[0] = 3.5;
> > > >
> > > > then the value is set. Otherwise it remains Neumann.
> > > >
> > > > - In the second case, one specifies that the boundary condition should
> > > > be u = g on the boundary for g a given Function. For a BDM element,
> > > > one would then say that the normal component of u should be the same
> > > > as the normal component of g on the boundary. Or do we want to say
> > > > that g should be the Dirichlet values of u, so that g would be a
> > > > scalar that returns the normal components in the case of BDM?
> > > >
> > >
> > > For more complicated elements, we could just impose Dirichlet boundary
> > > conditions weakly.
> > >
> > > Garth
> >
> > Agree.
> >
> > /Anders
> > _______________________________________________
> > DOLFIN-dev mailing list
> > DOLFIN-dev@xxxxxxxxxx
> > http://www.fenics.org/mailman/listinfo/dolfin-dev
> >
> 
>