On 24/05/17 12:12, Pan, Wei wrote:
> Hey guys,
>
>
> I'm having problems with the /Function.at()/ function. Below is an
> example code
>
>     from __future__ import division
>     from firedrake import *
>     import numpy as np
>
>     mesh_base = SquareMesh(8, 8, 1.)
>     mesh_hierarchy = MeshHierarchy(mesh_base, 6)
>
>     fine_msh = mesh_hierarchy[6]       /# 512 x 512 mesh/
>     fine_dg_space = FunctionSpace(fine_msh, "DG", 1)
>
>     x, y = SpatialCoordinate(fine_dg_space)
>     q0 = Function(fine_dg_space).interpolate(x)
>
>     x_array = np.array([_x / 512.0 for _x in range(513)])  # coordinates
>     # x_array[0] += 1. / 1000000.0 # a slight nudge away from boundary
>     also fails
>
>     for i in x_array:
>         for j in x_array:
>             print [i, j], q0.at(np.array([i, j]))
>
>
> My domain is [0, 1], mesh size 512x512. I keep on getting the
> following PointNotInDomainError when I run the above, where 0.43164062
> correspond to 221 of 512. What am I doing wrong?
>
>
>     [0.0, 0.427734375] 0.0
>     [0.0, 0.4296875] 0.0
>     [0.0, 0.431640625]
>     Traceback (most recent call last):
>       File "/home/wpan1/PycharmProjects/MyFirstProject/test_at.py",
>     line 22, in <module>
>         print [i, j], q0.at(np.array([i, j]))
>       File
>     "/home/wpan1/Documents/firedrake/src/firedrake/firedrake/function.py",
>     line 574, in at
>         raise PointNotInDomainError(self.function_space().mesh(),
>     points[i].reshape(-1))
>     firedrake.function.PointNotInDomainError: domain <Mesh #7> does
>     not contain point [ 0.          0.43164062]

The geometric predicates used to decide if a point is in a cell are
all done using floating point, so what might be happening is that this
point lies on the boundary between two cells, but to the specified
tolerance, it lies outside from both sides.

Does it help to relax the tolerance with which we check?

Say, for example:

q0.at(..., tolerance=1e-6)

Cheers,

Lawrence