Subject: [Ffc] Release plans

From: Anders Logg <[email-masked]>
Date: Fri, 29 Jan 2010 01:23:12 +0100

There seem to be just a couple of issues remaining in order of importance:

1. QuadratureElement

2. DOLFIN fem unit test

3. evaluate_basis_derivatives

4. RestrictedElement

Among these, I would say 1-2 are crucial to fix before 0.9.0,
but 3-4 are less crucial.

It would be good to have evaluate_basis_derivatives working, but I
don't know how much work is involved.

RestrictedElement is also good to have, but I would say it is less
important, especially since it's unclear what exactly it should do and
how it should be implemented.

We should also have the following in mind:

(i) The number of users relying on restricted elements is most likely
very small and they are very likely using the development version
anyway so if a fix for restricted elements comes next week, that will
be ok.

(ii) If we make a release of 0.9.0 now and package it for Debian/Ubuntu,
it is still possible to request Ubuntu to pull bug fixes for some
time. So we can release 0.9.0 now and then 0.9.1 in a couple of weeks.

One of the reasons I'm pressing on with the release is the interface
changes that we have made in both FFC and DOLFIN. It would be good to
introduce these as early as possible to avoid confusion about how to
write a form file (Coefficient, not Function) and how to implement
Expressions (Array, not std::vector).

--
Anders

---

From: "Garth N. Wells" <[email-masked]>
Date: Fri, 29 Jan 2010 00:26:32 +0000

> but 3-4 are less crucial.

Both 3 and 4 are crucial for me.

> Expressions (Array, not std::vector).

I think that it's all being rushed. With the major re-write of FFC it
should be given some time for user testing. I would be amazed if my
various solvers do not show up bugs.

Garth

> --
> Anders
> 
> 
> ------------------------------------------------------------------------
> 
> _______________________________________________
> Mailing list: https://launchpad.net/~ffc
> Post to     : [email-masked]
> Unsubscribe : https://launchpad.net/~ffc
> More help   : https://help.launchpad.net/ListHelp

---

From: Anders Logg <[email-masked]>
Date: Fri, 29 Jan 2010 09:32:02 +0100

> Both 3 and 4 are crucial for me.

I know, but I expect you are part of category (i) below, a limited
number of users who are running bzr tip anyway. :-)

> various solvers do not show up bugs.

Yes, it's being rushed, but there is time to release bug fixes.

To me, the minimum requirement is that we pass the regression tests in
FFC (perhaps except for evaluate_basis_derivatives and
ElementRestriction) and that we can compile and run all the DOLFIN
demos without problems. That should be enough for 0.9.0.

Even so, I'm not sure we will be able to manage to do that in time.
We should settle on something before lunch today. If we decide not to
press on with the release, we should take a good look at which
versions are currently packaged and see if those packages are in good
shape since those will go into the next Ubuntu LTS.

--
Anders

---

From: "Garth N. Wells" <[email-masked]>
Date: Fri, 29 Jan 2010 08:35:23 +0000

Anders Logg wrote:
>> Both 3 and 4 are crucial for me..
> 

Not exactly, because I share some solvers which require this
functionality with less adventurous collaborators who do use the packages.

Garth

---

From: Anders Logg <[email-masked]>
Date: Fri, 29 Jan 2010 09:54:37 +0100

> functionality with less adventurous collaborators who do use the packages.

If it's very important to you, could you help out with the design and
implementation of RestrictedElement. There seems to be some confusion
about what it should actually do. I still don't understand it.

--
Anders

---

From: "Garth N. Wells" <[email-masked]>
Date: Fri, 29 Jan 2010 09:31:54 +0000

Anders Logg wrote:
> implementation of RestrictedElement. 

I'd like to, but I won't be able to look at the code until next week.

> There seems to be some confusion

Here's an example of its use:

  http://www.dspace.cam.ac.uk/bitstream/1810/221687/1/solver.py

In the simplest case, it is just a restriction to cell facets. Kristian
would like to approach it more generally with restriction to cells,
measures an subdomains, but I think that there is a difference. What
RestrictedElement did is define functions that live only on the facets
of a cell (i.e. no internal dofs associated with basis functions that
are zero on facets. Restriction to measures are different, since the
function still lives on the entire cell, but happens to be evaluated on
a lower-dimensional entity, e.g. a surface passing through a cell.
Subdomains are different again and are just defining a subdomain on
which function is defined (a subset of cells).

Garth

---

From: Anders Logg <[email-masked]>
Date: Fri, 29 Jan 2010 10:57:37 +0100

> I'd like to, but I won't be able to look at the code until next week.

ok. Could Kristian comment on the possibility of finishing things up
today (without breaking yourself)? If it's not possible, it's not
possible and we should wait with the release 2-3 weeks so we have time
to test properly and fix a few more bugs. But perhaps it's possible to
fix almost everything, release 0.9.0 now and then 0.9.1 as a bug fix
release in a few weeks.

I think it's possible to ask Ubuntu to pull 0.9.1 as a bug fix for
0.9.0 in a few weeks, but it's likely they will not want to pull 0.9.0
as a major upgrade of 0.7.1.

>   http://www.dspace.cam.ac.uk/bitstream/1810/221687/1/solver.py

Looks good.

> which function is defined (a subset of cells).

Sub domains seem to be very different, but the other two cases just
seem to be a matter of some dofs being "active" and the other zeroed
out. This is what Marie suggested yesterday, that a restricted element
only considers a subset of the dofs of some given element.

The thing I don't understand yet is the selection of which dofs should
be active. If we think of the case with restriction to facets, then
the element needs to be restricted to different facets depending on
which facet we are integrating over, or are we always mapping one
specific facet of the reference cell to the current facet?

Say we have P1 elements in 2D which have 3 dofs. Then we could
restrict that element to the dofs on the first facet (facet 0). These
dofs are then labeled 1 and 2. But sometimes a facet in the mesh will
correspond to the edge between 0 and 1 or 0 and 2.

--
Anders

---

From: "Garth N. Wells" <[email-masked]>
Date: Fri, 29 Jan 2010 10:10:37 +0000

> only considers a subset of the dofs of some given element.

Sounds good.

> specific facet of the reference cell to the current facet?
It works the same way as the DG elements, just the internal dofs are
thrown away, which is the latter if the above, right?

> correspond to the edge between 0 and 1 or 0 and 2.

We don't restrict to individual facets, but to all facts of a cell.

Garth

---

From: Anders Logg <[email-masked]>
Date: Fri, 29 Jan 2010 11:28:07 +0100

> We don't restrict to individual facets, but to all facts of a cell.

That makes sense, but one thing still confuses me. Say that we have a
P1 element and restrict it to facets. Then all dofs are on the facets
so the result of the restriction is just a new P1 element. Same for P2
where the result again is a new P2 element. For P3, the result is P3
element minus just one dof. So does this make much difference for
other than very high degree elements?

--
Anders

---

From: "Garth N. Wells" <[email-masked]>
Date: Fri, 29 Jan 2010 10:31:05 +0000

Anders Logg wrote:
> where the result again is a new P2 element. 

Yes.

> For P3, the result is P3
> element minus just one dof.

Yes.

> So does this make much difference for

It is only needed for k > 2. It's important because everything in
FFC/UFL works for arbitrary orders.

Garth

---

From: Mehdi Nikbakht <[email-masked]>
Date: Fri, 29 Jan 2010 11:38:38 +0100

> only considers a subset of the dofs of some given element.

For the restriction with Measure(i.e. dc), it is not exactly the case.
The corresponding shape functions for thier dofs are defined in part of
cell and the rest are zeroed. We don't throw away any dofs in
compilation stage. 
Note that they might be inactive which will be determined in the solver.

For my case, I am using the restriction just as a marker to determine
discontinuous space.

Mehdi

---

From: Anders Logg <[email-masked]>
Date: Fri, 29 Jan 2010 11:51:37 +0100

> FFC/UFL works for arbitrary orders.

Yes, it should work for any order. I was just questioning the
usefulness of it if it results in standard P1 and P2 elements for
k = 1,2.

Then I think I understand how it all works.

But does the FFC demo make any sense? Could we simplify it so that it
just defines a finite element restricted to a facet and then some dS
integral?

It currently breaks the code we added for RestrictedElement yesterday
because it has a complex nesting of restriction, then mixed with
another element and then restricted again. It looks like this is not
how you are using it in your solver.

--
Anders

---

From: Anders Logg <[email-masked]>
Date: Fri, 29 Jan 2010 11:52:36 +0100

> For my case, I am using the restriction just as a marker to determine

Is that part currently broken and what needs to be added to make it
work?

--
Anders

---

From: "Garth N. Wells" <[email-masked]>
Date: Fri, 29 Jan 2010 10:54:54 +0000

Anders Logg wrote:
> But does the FFC demo make any sense? 

We have a demo?

> Could we simplify it so that it

Sounds good.

> how you are using it in your solver.

The concept is very simple, so its implementation should be pretty
simple too.

Garth

---

From: Anders Logg <[email-masked]>
Date: Fri, 29 Jan 2010 12:02:07 +0100

> > But does the FFC demo make any sense?

Yes, it's called ElementRestriction.ufl. :-)

> > just defines a finite element restricted to a facet and then some dS

Could you post something simple and I'll add it.

Yes.

--
Anders

---

From: "Garth N. Wells" <[email-masked]>
Date: Fri, 29 Jan 2010 12:00:24 +0000

> Yes, it's called ElementRestriction.ufl. :-)

I've just changed it. It doesn't test extensively, but it does hopefully
illustrate now better its purpose.

On a related note, it would be useful if we could have multiple forms in
one file and have all of them compiled. Can this be done already?

Garth

---

From: Kristian Oelgaard <[email-masked]>
Date: Fri, 29 Jan 2010 13:09:56 +0100 (CET)

> don't know how much work is involved.

First my apologies for the late reply, Internet was down again all morning and it pisses me off.

Now that BDM elements also compiles for gcc it looks like everything is in place for evaluate_basis.
CG and CR elements work, and so does all linear elements and DG constants. 2nd and 3rd order elements of BDM, RT, and DG type does not give the same result but there are explanations for this and at least it looks like they agree with the values you would get from FIAT if you tabulate on the reference element.

evaluate_basis_derivatives works for CG elements and for most elements of degree 1.

I'm working on some kind of unit testing framework for evaluate_basis and evaluate_basis_derivatives, such that we can make sure that at least cpp code and FIAT agrees on the reference element. If this test works, then the only bug can be in the transformation back to the physical element, but since this works for CG elements I don't see why it shouldn't for other elements too.
However, I'll put this development on hold until Quadrature and Restricted elements work.

> how it should be implemented.

It should be implemented like I did it in the old FFC :)
w.r.t. the demo, it should test if restrictions on a MixedElement are propagated to the sub elements and that a restricted element can appear in a MixedElement.

> time. So we can release 0.9.0 now and then 0.9.1 in a couple of weeks.

I think this sounds like a good plan. I'll work on getting QuadratureElement up and running and the RestrictedElement.
If there are any bugs in evaluate_basis_derivatives we'll release it in 0.9.1.
The code really needs some cleaning up, and we should also rearrange some of the utility functions, but that can wait until 0.9.1.

Kristian

> -----BEGIN PGP SIGNATURE-----
> Version: GnuPG v1.4.9 (GNU/Linux)
>
> iEYEARECAAYFAktiKnAACgkQTuwUCDsYZdHEmgCcD3jSWQfDtVnNktS2+XMgQMhg
> FkQAoJL6rELvZTzC+8fea5uMIYohDwUV
> =aMac
> -----END PGP SIGNATURE-----
>

---

From: Marie Rognes <[email-masked]>
Date: Fri, 29 Jan 2010 13:11:36 +0100

Kristian Oelgaard wrote:
> CG and CR elements work, and so does all linear elements and DG 
> constants. 2nd and 3rd order elements of BDM, RT, and DG type does not 
> give the same result but there are explanations for this and at least 
> it looks like they agree with the values you would get from FIAT if 
> you tabulate on the reference element.
>
> I'm working on some kind of unit testing framework for evaluate_basis 
> and evaluate_basis_derivatives, such that we can make sure that at 
> least cpp code and FIAT agrees on the reference element. If this test 
> works, then the only bug can be in the transformation back to the 
> physical element, but since this works for CG elements I don't see why 
> it shouldn't for other elements too.
> w.r.t. the demo, it should test if restrictions on a MixedElement are 
> propagated to the sub elements and that a restricted element can 
> appear in a MixedElement.
>

I'm taking care of RestrictedElement (onto facets/cells etc) at ze moment.

--
Marie

> I think this sounds like a good plan. I'll work on getting 
> QuadratureElement up and running and the RestrictedElement.
> The code really needs some cleaning up, and we should also rearrange 
> some of the utility functions, but that can wait until 0.9.1.
>
> Kristian
>

---

From: Kristian Oelgaard <[email-masked]>
Date: Fri, 29 Jan 2010 13:13:32 +0100 (CET)

On 29 January 2010 13:11, Marie Rognes <[email-masked]> wrote:
>> CG and CR elements work, and so does all linear elements and DG constants.
>> 2nd and 3rd order elements of BDM, RT, and DG type does not give the same
>> result but there are explanations for this and at least it looks like they
>> agree with the values you would get from FIAT if you tabulate on the
>> reference element.
>>
>> I'm working on some kind of unit testing framework for evaluate_basis and
>> evaluate_basis_derivatives, such that we can make sure that at least cpp
>> code and FIAT agrees on the reference element. If this test works, then the
>> only bug can be in the transformation back to the physical element, but
>> since this works for CG elements I don't see why it shouldn't for other
>> elements too.
>> propagated to the sub elements and that a restricted element can appear in a
>> MixedElement.
>>
>

OK, I'll do Quadrature then.

Kristian

>> The code really needs some cleaning up, and we should also rearrange some
>> of the utility functions, but that can wait until 0.9.1.
>>
>> Kristian
>>

---

From: Marie Rognes <[email-masked]>
Date: Fri, 29 Jan 2010 16:55:19 +0100

Marie Rognes wrote:
>> constants. 2nd and 3rd order elements of BDM, RT, and DG type does 
>> not give the same result but there are explanations for this and at 
>> least it looks like they agree with the values you would get from 
>> FIAT if you tabulate on the reference element.
>>
>> physical element, but since this works for CG elements I don't see 
>> why it shouldn't for other elements too.

This was actually slightly more mind-boggling than planned.
I ended up letting FFC convert restricted mixed elements to mixed 
restricted elements.
That kept the codegeneration complexity low and my head from exploding.

I'll clean up and push tomorrow, unless that is too late ;)

Also, question for Kristian: About the IR for evaluate_basis for mixed 
elements: You want a (non-nested) list of sub-elements, right?

--
Marie

---

From: Kristian Oelgaard <[email-masked]>
Date: Fri, 29 Jan 2010 17:08:58 +0100 (CET)

On 29 January 2010 16:55, Marie Rognes <[email-masked]> wrote:
>>> constants. 2nd and 3rd order elements of BDM, RT, and DG type does not give
>>> the same result but there are explanations for this and at least it looks
>>> like they agree with the values you would get from FIAT if you tabulate on
>>> the reference element.
>>>
> I ended up letting FFC convert restricted mixed elements to mixed restricted
> elements.

Yes, I think that should work. But why do you need to bother if you create a 'normal' MixedElement? Then nothing should change with respect to the functions in representation.py

Kristian

---

From: Marie Rognes <[email-masked]>
Date: Fri, 29 Jan 2010 18:46:48 +0100

Kristian Oelgaard wrote:
>>>> Now that BDM elements also compiles for gcc it looks like 
>>>> everything is
>>>> in place for evaluate_basis.
>>>> not give
>>>> evaluate_basis_derivatives works for CG elements and for most 
>>>> elements of
>>>> degree 1.
>>>>
>>>> appear in a
>> restricted
> Yes, I think that should work. But why do you need to bother if you 
> create a 'normal' MixedElement? Then nothing should change with 
> respect to the functions in representation.py
>

The IR for evaluate_basis for mixed element iterate[s|d] over the 
ufl_element. I don't want to do that (since that would be an 
ufl.ElementRestriction.)

--
Marie

---

From: Kristian Oelgaard <[email-masked]>
Date: Fri, 29 Jan 2010 19:01:27 +0100 (CET)

On 29 January 2010 18:46, Marie Rognes <[email-masked]> wrote:
>>>>> give
>>>>> evaluate_basis_derivatives works for CG elements and for most elements
>>>>> of
>>>>> degree 1.
>>>>>
>>>>> and
>> Yes, I think that should work. But why do you need to bother if you create
>> a 'normal' MixedElement? Then nothing should change with respect to the
>> functions in representation.py
>>
>

I don't follow. If you have a mixed element

mixed_element = elem0 + elem1

and then restrict it to facets:

element = mixed_element["facet"]

which is now a RestrictedElement. Then, as I understand it, you create a new MixedElement equivalent to:

element = elem0["facet"] + eleme1["facet"]

in this case the evaluate_basis IR (and all other IRs) will only see a 'normal' MixedElement and we can proceed without problems.

Kristian

---

From: Marie Rognes <[email-masked]>
Date: Fri, 29 Jan 2010 22:02:22 +0100

Kristian Oelgaard wrote:
>>>>>> like they agree with the values you would get from FIAT if you 
>>>>>> tabulate
>>>>>> on
>>>>>> the reference element.
>>>>>>
>>>>>> elements
>>>>>> code and FIAT agrees on the reference element. If this test 
>>>>>> works, then
>>>>>> the
>>>>>> appear
>>>>>> in a
>>> create
> which is now a RestrictedElement. 

on the UFL side, yes.

> Then, as I understand it, you create a new MixedElement equivalent to:
>

on the FFC side, yes

> in this case the evaluate_basis IR (and all other IRs) will only see a 
> 'normal' MixedElement and we can proceed without problems.
>

Exactly. That is the bonus of this approach :)

But I did not want to tamper with the UFL representation of the element. 
Hence, it is way easier to iterate
over the FFC element rather than the UFL element in the IR of 
evaluate_basis. Just wanted
to check with you that the way I planned doing that was correct. For the 
other IRs I know what
the correct way of doing this is ;)

--
Marie