dolfin team mailing list archive

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

On 23 November 2011 09:37, Kristian Ølgaard <k.b.oelgaard@xxxxxxxxx> wrote:
> On 22 November 2011 21:15, Kristian Ølgaard <k.b.oelgaard@xxxxxxxxx> wrote:
>> On 22 November 2011 21:07, Anders Logg <logg@xxxxxxxxx> wrote:
>>> On Tue, Nov 22, 2011 at 07:34:38PM +0000, Garth N. Wells wrote:
>>>> On 22 November 2011 19:12, Kristian Ølgaard <k.b.oelgaard@xxxxxxxxx> wrote:
>>>> > Hi,
>>>> >
>>>> > The IntervalCell::facet_area in DOLFIN currently returns 0.0, while
>>>> > the facet determinant (a scale factor) 'det = 1.0;'
>>>> > in the generated code for the tabulate_tensor function in 1D (facet integrals).
>>>> >
>>>> > Is 0.0 the correct return value, or should it return 1.0?
>>>> >
>>>>
>>>> I would think 1.0.
>>
>> That's what I'm leaning towards too.
>>
>>> Facet is the thing of dimension one lower than the cell dimension so
>>> in this case 1 - 1 so it would mean the area of an endpoint. So I
>>> believe 0 is correct.
>>
>> I see your 'point' :), but is it consistent with how we evaluate
>> integrals on endpoints?
>
> Running this script:
>
> from dolfin import *
>
> cells = [interval, triangle, tetrahedron]
> meshes = [UnitInterval(1), UnitSquare(1,1), UnitCube(1,1,1)]
>
> for cell, mesh in zip(cells, meshes):

btw, you can do
    cell = mesh.ufl_cell()

>    c = Constant(1, cell)
>    print
>    print cell
>    print "volume:  ", assemble(c*dx, mesh=mesh)
>    print "surface: ", assemble(c*ds, mesh=mesh)
>
> produces:
>
> <interval cell in R1>
> volume:   1.0
> surface:  2.0
>
> <triangle cell in R2>
> volume:   1.0
> surface:  4.0
>
> <tetrahedron cell in R3>
> volume:   1.0
> surface:  6.0
>
> If we want surface = 0.0 for the interval, 'det' must be set to zero
> in the generated code.
> However, doing so will mean that f*ds etc. becomes zero too, which is
> not very useful.
>
> Kristian

As Kristian shows here, I think defining it to 1.0 seems more useful,
and consistent with the definition of 1.0*ds.

Martin