← Back to team overview

ffc team mailing list archive

Re: "Bug" in evaluate_basis

 



Anders Logg wrote:
@Garth: Yes, they are defined everywhere. For example, v1 = 1 - x - y,
v2 = x and v3 = y are a nodal basis for the reference triangle, but
they are also a basis for P1 on R^2. You can evaluate 1 - x - y for
any values of x and y.


I don't agree at all. The function 1-x-y can be evaluated anywhere, but the basis functions are defined only on the cell. If the basis functions were defined everywhere, the method would loose all sparsity.

Garth

@Matt: We use the local basis on a triangle K as a basis for P2 (and
other spaces) on a patch of elements surrounding K. We make an Ansatz
in terms of the local basis functions. So it's not for extrapolation
of a computed solution.

We could also use the Bernstein polynomials on K or any other
reasonable well-conditioned basis, but it would be convenient to use
the basis that's already in there.

--
Anders


On Wed, Sep 16, 2009 at 09:05:20AM -0500, Matthew Knepley wrote:
Isn't this strange since the error bounds no longer apply outside the element
(it is no longer interpolation, but extrapolation)?

  Matt

On Wed, Sep 16, 2009 at 8:59 AM, Anders Logg <logg@xxxxxxxxx> wrote:

    Marie and I came across a "bug" when trying to use evaluate_basis in a
    way it was probably not intended for.

    In particular, we want to evaluate the basis functions for a cell in a
    point that lies outside the cell. This fails when we evaluate the
    basis functions in a point where the collapsing square map in FIAT is
    singular.

    Here's a simple script that illustrates the problem:

     from dolfin import *
     import numpy

     mesh = UnitSquare(1, 1, "left")
     V = FunctionSpace(mesh, "CG", 1)
     element = V.dolfin_element()
     cell = Cell(mesh, 0)
     value = numpy.zeros(1)

     # Value should be 1.0 (constant along x = 1)
     element.evaluate_basis(1, value, numpy.array((1.0, 1.0000000000001)),
    cell)
     print value

     # Should also be 1.0 but result is 0 (eps)
     element.evaluate_basis(1, value, numpy.array((1.0, 1.0000000000000)),
    cell)
     print value

    The problem is the following lines in evaluate_basis:

     // Map coordinates to the reference square
     if (std::abs(y - 1.0) < 1e-14)
       x = -1.0;
     else
       x = 2.0 *x/(1.0 - y) - 1.0;
     y = 2.0*y - 1.0;

    In the above example, y = 1 so x is mapped to -1.0 on Rob's reference
    square.

    Is there a simple fix that would allow us to evaluate the basis
    functions outside the element?


    -----BEGIN PGP SIGNATURE-----
    Version: GnuPG v1.4.9 (GNU/Linux)

    iEYEARECAAYFAkqw71UACgkQTuwUCDsYZdEIagCeK7i8b18hMN/aD1tQioosHu6L
    GVsAn28SJ86qU+zUm2AUuYLpvTgTD7ds
    =qOFA
    -----END PGP SIGNATURE-----

    _______________________________________________
    FFC-dev mailing list
    FFC-dev@xxxxxxxxxx
    http://www.fenics.org/mailman/listinfo/ffc-dev






_______________________________________________
FFC-dev mailing list
FFC-dev@xxxxxxxxxx
http://www.fenics.org/mailman/listinfo/ffc-dev


------------------------------------------------------------------------

_______________________________________________
FFC-dev mailing list
FFC-dev@xxxxxxxxxx
http://www.fenics.org/mailman/listinfo/ffc-dev




Follow ups

References