dolfin team mailing list archive

[Kristian Ølgaard <k.b.oelgaard@xxxxxxxxx>]

On 25 November 2011 12:56, Anders Logg <logg@xxxxxxxxx> wrote:
> Looks like you need to check that you are not running in parallel in
> that test since DOLFIN will otherwise try to partition the mesh which
> doesn't work if you only have one cell.

Done.

Kristian

> --
> Anders
>
>
> On Fri, Nov 25, 2011 at 10:58:17AM +0100, Kristian Ølgaard wrote:
>> On 24 November 2011 20:33, Anders Logg <logg@xxxxxxxxx> wrote:
>> > On Thu, Nov 24, 2011 at 01:43:20PM +0100, Kristian Ølgaard wrote:
>> >> On 23 November 2011 10:54, Martin Sandve Alnæs <martinal@xxxxxxxxx> wrote:
>> >> > 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.
>> >>
>> >> I've pushed the change to DOLFIN, the output of the extended script:
>> >>
>> >> from dolfin import *
>> >>
>> >> for mesh in [UnitInterval(1), UnitSquare(1,1), UnitCube(1,1,1)]:
>> >>     c = Constant(1, mesh.ufl_cell())
>> >>     c0 = mesh.ufl_cell().facet_area
>> >>     c1 = FacetArea(mesh)
>> >>     print
>> >>     print mesh.ufl_cell()
>> >>     print "volume:  ", assemble(c*dx, mesh=mesh)
>> >>     print "surface: ", assemble(c*ds, mesh=mesh)
>> >>     print "facet_area ufl/ffc: ", assemble(c0*ds, mesh=mesh)
>> >>     print "FacetArea  dolfin:  ", assemble(c1*ds, mesh=mesh)
>> >>
>> >> is
>> >>
>> >> <interval cell in R1>
>> >> volume:   1.0
>> >> surface:  2.0
>> >> facet_area ufl/ffc:  2.0
>> >> FacetArea  dolfin:   2.0
>> >>
>> >> <triangle cell in R2>
>> >> volume:   1.0
>> >> surface:  4.0
>> >> facet_area ufl/ffc:  4.0
>> >> FacetArea  dolfin:   4.0
>> >>
>> >> <tetrahedron cell in R3>
>> >> volume:   1.0
>> >> surface:  6.0
>> >> facet_area ufl/ffc:  3.0
>> >> FacetArea  dolfin:   3.0
>> >>
>> >> Kristian
>> >
>> > ok.
>> >
>> > Can you make it into a unit test? (Maybe you already did, I haven't checked.)
>>
>> I put the above script in
>> test/unit/function/python/SpecialFunctions.py (not pushed to repo)
>>
>> but doing
>>
>> mpirun -np 3  python ./SpecialFunctions.py
>>
>> results in the error:
>>
>> Process 0: Number of global vertices: 2
>> Process 0: Number of global cells: 1
>> EE
>> ======================================================================
>> ERROR: testFacetArea (__main__.SpecialFunctions)
>> ----------------------------------------------------------------------
>> Traceback (most recent call last):
>>   File "./SpecialFunctions.py", line 29, in testFacetArea
>>     references = [(UnitInterval(1), 2, 2),]#\
>>   File "/home/oelgaard/software/fenics/lib/python2.7/site-packages/dolfin/cpp.py",
>> line 12959, in __init__
>>     _cpp.UnitInterval_swiginit(self,_cpp.new_UnitInterval(*args))
>> RuntimeError:
>>
>> *** -------------------------------------------------------------------------
>> *** DOLFIN encountered an error. If you are not able to resolve this issue
>> *** using the information listed below, you can ask for help at
>> ***
>> ***     https://answers.launchpad.net/dolfin
>> ***
>> *** Remember to include the error message listed below and, if possible,
>> *** include a *minimal* running example to reproduce the error.
>> ***
>> *** -------------------------------------------------------------------------
>> *** Error:  Unable to compute mesh partitioning using ParMETIS.
>> *** Reason: ParMETIS cannot be used if a process has no cells (graph
>> nodes). Use SCOTCH to perform partitioning instead.
>> *** Where:  This error was encountered inside ParMETIS.cpp.
>> *** -------------------------------------------------------------------------
>>
>>
>> ----------------------------------------------------------------------
>> Ran 1 test in 0.003s
>>
>> FAILED (errors=1)
>>
>> ======================================================================
>> ERROR: testFacetArea (__main__.SpecialFunctions)
>> ----------------------------------------------------------------------
>> Traceback (most recent call last):
>>   File "./SpecialFunctions.py", line 29, in testFacetArea
>>     references = [(UnitInterval(1), 2, 2),]#\
>>   File "/home/oelgaard/software/fenics/lib/python2.7/site-packages/dolfin/cpp.py",
>> line 12959, in __init__
>>     _cpp.UnitInterval_swiginit(self,_cpp.new_UnitInterval(*args))
>> RuntimeError:
>>
>> *** -------------------------------------------------------------------------
>> *** DOLFIN encountered an error. If you are not able to resolve this issue
>> *** using the information listed below, you can ask for help at
>> ***
>> ***     https://answers.launchpad.net/dolfin
>> ***
>> *** Remember to include the error message listed below and, if possible,
>> *** include a *minimal* running example to reproduce the error.
>> ***
>> *** -------------------------------------------------------------------------
>> *** Error:  Unable to compute mesh partitioning using ParMETIS.
>> *** Reason: ParMETIS cannot be used if a process has no cells (graph
>> nodes). Use SCOTCH to perform partitioning instead.
>> *** Where:  This error was encountered inside ParMETIS.cpp.
>> *** -------------------------------------------------------------------------
>>
>>
>> ----------------------------------------------------------------------
>> Ran 1 test in 0.003s
>>
>> FAILED (errors=1)
>>
>>
>> after which the process hangs.
>>
>> Kristian
>>
>> >
>