fiat team mailing list archive
-
fiat team
-
Mailing list archive
-
Message #00303
Re: [FFC-dev] "Bug" in evaluate_basis
@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.
@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
Attachment:
signature.asc
Description: Digital signature
Follow ups
References