Issue #10 resolved

Erroneous output when tdim != gdim for compute_vertex_values on RT

Chris Richardson
created an issue

See https://bitbucket.org/fenics-project/dolfin/issue/75

FFC generates code for compute_vertex_values() which is the wrong dimension when tdim != gdim (at least for RT elements). See example in dolfin issue #75.

Comments (12)

  1. Martin Alnæs

    I think the bug is somewhere in this behaviour in fiatinterface.py:

    def reference_cell(dim):
        ...
        return FIAT.ufc_simplex(cellname2dim[dim])
    
    # These lines in _create_fiat_element:
    fiat_cell = reference_cell(cell.cellname())
    element = ElementClass(fiat_cell, degree)
    

    Note how gdim is lost through cellname2dim[cell.cellname()] and FIAT builds a pure 2D RT element. Now the value size of element is 2 instead of 3 which results in the bug you tried to fix with your patch. I'm not sure what the proper solution is.

  2. Martin Alnæs

    For the record, there are probably more issues with this code for uncommon combinations such as mixed and fancy elements on manifolds. I added a safeguard in this particular code that should detect if the same problem pops up with another combo.

  3. Martin Alnæs

    I had to make another fix to pass the tests: when the safeguard fails, interpolate_vertex_values is generated to throw an exception instead of failing at runtime. This will happen e.g. with a RT1*DG0 mixed element on a manifold (see the interpolate_vertex_values implementations in ProjectionManifold.h). Fixing this looks a bit trickier (because of code organization), I have no plans to do so.

  4. Marie Rognes

    The original issue (DOLFIN issue 75, which was marked resolved, but was definitely still failing) now has a fix in DOLFIN next. I'll resolve this when that fix has graduated to master.

  5. Marie Rognes

    For completeness here: Martin Alnæs, FIAT should build an element of dimension tdim (topological dimension) as the reference element does not need to know about the geometrical dimension of the physical element. I'm not sure what the fixes made in the fall were (and I'm not sure if I want to ...)

  6. Log in to comment