Assembly on manifolds

On Thu, Sep 27, 2012 at 12:01:32PM -0000, Charl wrote:
> New question #209712 on DOLFIN:
> https://answers.launchpad.net/dolfin/+question/209712
>
> Hello,
>
> Just a quick question about the current status of this feature. The FEniCS book mentions the advantages that this would have, e.g. support for assembly over R^2 elements (triangles) which are embedded in R^3.
>
> There are several questions relating to this, e.g.
> https://answers.launchpad.net/dolfin/+question/193816
> https://answers.launchpad.net/dolfin/+question/169094
> ... but the youngest is dated beginning of this year.
>
> Has there been any progress in implementing this feature? Is it part
> of a development version that one can download and experiment with?
> Otherwise, is it possible to 'hack' this together in some other way?

No progress that I know of but there is interest in adding support for
this so it might come (can't promise when). So the best option now
would be to hack something together.

--
Anders

On 09/27/2012 02:11 PM, Anders Logg wrote:
> Question #209712 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/209712
>
> Status: Open => Answered
>
> Anders Logg proposed the following answer:
> On Thu, Sep 27, 2012 at 12:01:32PM -0000, Charl wrote:
>> New question #209712 on DOLFIN:
>> https://answers.launchpad.net/dolfin/+question/209712
>>
>> Hello,
>>
>> Just a quick question about the current status of this feature. The FEniCS book mentions the advantages that this would have, e.g. support for assembly over R^2 elements (triangles) which are embedded in R^3.
>>
>> There are several questions relating to this, e.g.
>> https://answers.launchpad.net/dolfin/+question/193816
>> https://answers.launchpad.net/dolfin/+question/169094
>> ... but the youngest is dated beginning of this year.
>>
>> Has there been any progress in implementing this feature? Is it part
>> of a development version that one can download and experiment with?
>> Otherwise, is it possible to 'hack' this together in some other way?
> No progress that I know of but there is interest in adding support for
> this so it might come (can't promise when). So the best option now
> would be to hack something together.
>

I in particular have definite interest, and plan on working on this
in the course of the (late) fall.

--
Marie

--
Marie Elisabeth Rognes
Senior Research Scientist
Simula Research Laboratory
P.O. Box 134, 1325 Lysaker
http://home.simula.no/~meg/

Hello,

Great to see the shared interest. I would however like to start tinkering with this sometime soon.

Would an option be to take the fully assembled system matrix, and force all elements to zero that do not correspond to facets within the 2D manifold? This would essentially sever the manifold vertices from the surrounding vertices. I'm not sure however this would be numerically correct. (And wasted computation set aside.)

Any alternative suggestions are warmly appreciated!

It rather depends what sort of equation you are solving and therefore what sorts of function spaces you are dealing with. If you only have scalar fields then something like what you propose or the following might work. If you have vector fields then unless you have proper manifold fields you will run into trouble with keeping the vectors tangent to the manifold.

If you can arrange for the manifold you are interested in to be the surface of a 3D mesh (or part thereof) then you might be able to assemble your equations by simply composing forms integrated over ds. I think that probably gets you the right matrix entries. You'd then have to do some matrix surgery to deal with all the zero rows. This probably involves diving straight into the matrix backend to (for example) replace zeros on the matrix diagonal with 1. Not very memory efficient but it shouldn't change the convergence rate of an iterative solver.

Because there is no measure which would enable you to integrate along an edge of a 3D mesh (and there is also currently no infrastructure for getting the normals right), there isn't an easy way to do integrals over the facets of the manifold elements (ie over the edges of the triangles).

On Thu, Sep 27, 2012 at 03:56:05PM -0000, David Ham wrote:
> Question #209712 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/209712
>
> David Ham posted a new comment:
> It rather depends what sort of equation you are solving and therefore
> what sorts of function spaces you are dealing with. If you only have
> scalar fields then something like what you propose or the following
> might work. If you have vector fields then unless you have proper
> manifold fields you will run into trouble with keeping the vectors
> tangent to the manifold.
>
> If you can arrange for the manifold you are interested in to be the
> surface of a 3D mesh (or part thereof) then you might be able to
> assemble your equations by simply composing forms integrated over ds. I
> think that probably gets you the right matrix entries. You'd then have
> to do some matrix surgery to deal with all the zero rows. This probably
> involves diving straight into the matrix backend to (for example)
> replace zeros on the matrix diagonal with 1. Not very memory efficient
> but it shouldn't change the convergence rate of an iterative solver.

The matrix surgery can be done with A.ident_zeros() but it would be
inefficient. Something else that might help for this is upcoming
support for restricted function spaces (so you can restrict the 3D
function space to the surface), which I plan to work on in the near
future.

--
Anders

> Because there is no measure which would enable you to integrate along an
> edge of a 3D mesh (and there is also currently no infrastructure for
> getting the normals right), there isn't an easy way to do integrals over
> the facets of the manifold elements (ie over the edges of the
> triangles).
>