ufl team mailing list archive

["Martin Sandve Alnæs" <martinal@xxxxxxxxx>]

2008/6/19 Kent-Andre Mardal <kent-and@xxxxxxxxx>:
>
> tor, 19.06.2008 kl. 12.17 +0200, skrev Martin Sandve Alnæs:
>> 2008/6/19 Kent-Andre Mardal <kent-and@xxxxxxxxx>:
>> >
>> > tor, 19.06.2008 kl. 11.38 +0200, skrev Martin Sandve Alnæs:
>> >> 2008/6/19 Kent-Andre Mardal <kent-and@xxxxxxxxx>:
>> >> >
>> >> > Should we include Dirichlet conditions in UFL ?
>> >> >
>> >> > It could be done with the notion of traces:
>> >> >
>> >> > g = Function(...)
>> >> > u = TrialFunction(...)
>> >> > T = FacetTrace(...)
>> >>
>> >> From a programming point of view, this line makes no sense to me:
>> >>
>> >> > T*u = g
>> >>
>> >> What it means mathematically is clear, but this can never be valid
>> Python.
>> >
>> > You're right. Maybe
>> > T*u == g.
>> > Is this possible in Python ? It is in C++...
>>
>> The assignment operator doesn't exist in Python.
>>
>> Even if it did work, say in C++ or using ==, "T*u" would a a temporary
>> object,
>> and "(T*u) == g" would be a call to (T*u).operator==(g) returning
>> another
>> temporary object, which would then dissapear on the next line.
>>
>>
>> >> Do you mean something like:
>> >>
>> >> a = u*v*dx + dot(u,n)*v*ds(0) + g*T(1)
>> >>
>> >> where T(i) and ds(i) would refer to the same boundary?
>> >>
>> >
>> > Could be but I think this is less mathematically clear
>> > (for one thing a would not be bilinear)
>> >
>> > c = T*u == g.
>> > a = u*v*dx + dot(u,n)*v*ds(0)
>>
>> So you don't want the BC to be part of the form.
>
> It is a form, just not a bilinear one.

UFL will have, as FFC has, operators "lhs(form)" and "rhs(form)"
which extracts the bilinear and linear integrals from a form.
It is possible that these could "do the right thing" with a form
like the one above including DBC.

We could collect them more explicitly:
  a = u*v*dx
  L = f*v*ds
  bc = T*u == g
  Axb = System(a, L, bc)


>> Where should it be used? Why should it be part of UFL
>> and not just a dolfin component?
>>
>
> Could be that we should leave Dirichlet bc as it is. It's just
> that T is in some respect very similar to ds.

Don't misunderstand me, I think it is a good idea,
I'm just trying to filter out all problems I see.

We must keep in mind that dolfin is a C++ library,
UFL is a Python implementation, and dolfin should
only relate to UFC, not UFL directly. The exception
is perhaps a thin layer in PyDOLFIN where jit is called
from a form compiler to convert from UFL to UFC.
It is important to stick with the interfaces we agree on.


>> If we are to include this in UFL, it should support arbitrary UFL expressions
>> in place of g, and we would also have to support component-wise conditions:
>>
>> just_some_expression = (dot(g, n) + inner(grad(g), grad(g)))
>> c = T*(u[0]) == just_some_expression
>
> This makes sense, eg in elasticity you have Neumann in one component and
> Dirichlet in another on the same sub-domain.
>
>>
>> This means we must generate code for it, and need to extend UFC to include it.
>>
>> If we only support the basic "u = g" condition, I don't see the point
>> with adding it to UFL.
>>
>
> I am not sure how much this affect UFC, but it could for instance be
> nice to write
>
> T*u*normal = g
>
> where g is a Function and normal is the outward normal.

(As above this syntax isn't possible in the embedded language approach
we're using.
If we designed our own language from scratch (like mython?), that
could be possible.
Basically, we need to embed the information we wish to pass around into objects,
restricted by the syntax Python provides.)

That aside, we would have to define something flexible enough, which
may be difficult.
Basically, we could do something like:

  bc = T*u_components == expression

Here expression is any UFL-expression, which is compiled into some
future extension of UFC.

u_components will typically be:
  u - all values (in case of non-scalar element)
  u[0] - a single component (vector valued element, must be fixed integer)
  u[0, 1] - a single component (tensor valued element, must be fixed integers)
  dot(u, n) - normal component of u (handled differently depending on element!)

We _could_ define a small set of valid expressions for u_components
like the list above, if we are sure we can cover everything we need.

The extensions we need in UFC will be something to evaluate the
expression above given coefficients (similar to the integral classes
and tabulate_tensor), and something to tabulate which dofs to touch.

How to handle the dot(u, n), I'm not sure?


Regarding the subdomain numbering: we may want to have separate
numbering of T(i) and ds(i), to allow for BCs with some components
Dirichlet and some Neumann? Any examples of this?

Also, in many cases you'll want to just fix one point to avoid translation.
This kind of situation makes it difficult...


> But anyway, this is just a thought that has not been thought properly
> through yet...

Those are the best kind!

--
Martin