Dear all, When trying to evaluate a function u on some mesh points with mesh = Mesh('coarse/Omega_coarse.msh') V = FunctionSpace(mesh, "CG",1) V_Vec = VectorFunctionSpace(mesh, "CG", 1) u = Function(V) points = Function(V_vec).interpolate(SpatialCoordinate(mesh)) dat = u.at(points.dat.data_ro) I get the following error: File "lin_dyn_unicell.py", line 117, in <module> dat = u.at(points.dat.data_ro) File "/Users/gidel/Documents/Firedrake/firedrake/src/firedrake/firedrake/function.py", line 585, in at raise PointNotInDomainError(self.function_space().mesh(), points[i].reshape(-1)) firedrake.function.PointNotInDomainError: domain <Mesh #1> does not contain point [200. 66.66666667] How is that possible while both function spaces are defined on the same mesh? Note that in this example I am using the same domain to make sure the points are the same in the two cases, but I used to run this with different meshes and it worked on my laptop with an old firedrake version but not on my new machine on which I installed Firedrake yesterday. Thank you for your feedback, Floriane
Dear Floriane,
On 26 Oct 2018, at 11:06, Floriane Gidel [RPG] <mmfg@leeds.ac.uk> wrote:
Dear all,
When trying to evaluate a function u on some mesh points with
mesh = Mesh('coarse/Omega_coarse.msh') V = FunctionSpace(mesh, "CG",1) V_Vec = VectorFunctionSpace(mesh, "CG", 1) u = Function(V) points = Function(V_vec).interpolate(SpatialCoordinate(mesh)) dat = u.at(points.dat.data_ro)
I get the following error:
File "lin_dyn_unicell.py", line 117, in <module> dat = u.at(points.dat.data_ro) File "/Users/gidel/Documents/Firedrake/firedrake/src/firedrake/firedrake/function.py", line 585, in at raise PointNotInDomainError(self.function_space().mesh(), points[i].reshape(-1)) firedrake.function.PointNotInDomainError: domain <Mesh #1> does not contain point [200. 66.66666667]
How is that possible while both function spaces are defined on the same mesh?
You have run into the delightful world of computational geometry. Unfortunately points might not be determined to be in any cell if they fall close to cell boundaries. You can reduce the tightness of the tolerance used to check if a point is in a cell by saying: data = u.at(points, tolerance=1e-6) Thanks, Lawrence
Ok thanks Lawrence. ________________________________ De : firedrake-bounces@imperial.ac.uk <firedrake-bounces@imperial.ac.uk> de la part de Lawrence Mitchell <wence@gmx.li> Envoyé : vendredi 26 octobre 2018 11:12 À : firedrake Objet : Re: [firedrake] .at() and PointNotInDomainError Dear Floriane,
On 26 Oct 2018, at 11:06, Floriane Gidel [RPG] <mmfg@leeds.ac.uk> wrote:
Dear all,
When trying to evaluate a function u on some mesh points with
mesh = Mesh('coarse/Omega_coarse.msh') V = FunctionSpace(mesh, "CG",1) V_Vec = VectorFunctionSpace(mesh, "CG", 1) u = Function(V) points = Function(V_vec).interpolate(SpatialCoordinate(mesh)) dat = u.at(points.dat.data_ro)
I get the following error:
File "lin_dyn_unicell.py", line 117, in <module> dat = u.at(points.dat.data_ro) File "/Users/gidel/Documents/Firedrake/firedrake/src/firedrake/firedrake/function.py", line 585, in at raise PointNotInDomainError(self.function_space().mesh(), points[i].reshape(-1)) firedrake.function.PointNotInDomainError: domain <Mesh #1> does not contain point [200. 66.66666667]
How is that possible while both function spaces are defined on the same mesh?
You have run into the delightful world of computational geometry. Unfortunately points might not be determined to be in any cell if they fall close to cell boundaries. You can reduce the tightness of the tolerance used to check if a point is in a cell by saying: data = u.at(points, tolerance=1e-6) Thanks, Lawrence _______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
It should also be pointed out that you might (depending on circumstances) be making life overly complicated for yourself. The mesh coordinates are just a Function themselves, so if you have a CG1 Function then u.dat.data_ro is exactly the Function evaluated at the mesh vertices, in the same order as mesh.coordinate.dat.data_ro . Furthermore, you can achieve the same thing for the nodes of any CG or DG space by reinterpolating the coordinates. See the instructions at: https://www.firedrakeproject.org/interpolation.html#interpolation-from-exter... (admittedly those instructions are for interpolating into Firedrake rather than out of, but the same principle applies) Regards, David From: <firedrake-bounces@imperial.ac.uk> on behalf of "Floriane Gidel [RPG]" <mmfg@leeds.ac.uk> Reply-To: firedrake <firedrake@imperial.ac.uk> Date: Friday, 26 October 2018 at 13:07 To: firedrake <firedrake@imperial.ac.uk> Subject: Re: [firedrake] .at() and PointNotInDomainError Ok thanks Lawrence. ________________________________ De : firedrake-bounces@imperial.ac.uk <firedrake-bounces@imperial.ac.uk> de la part de Lawrence Mitchell <wence@gmx.li> Envoyé : vendredi 26 octobre 2018 11:12 À : firedrake Objet : Re: [firedrake] .at() and PointNotInDomainError Dear Floriane,
On 26 Oct 2018, at 11:06, Floriane Gidel [RPG] <mmfg@leeds.ac.uk> wrote:
Dear all,
When trying to evaluate a function u on some mesh points with
mesh = Mesh('coarse/Omega_coarse.msh') V = FunctionSpace(mesh, "CG",1) V_Vec = VectorFunctionSpace(mesh, "CG", 1) u = Function(V) points = Function(V_vec).interpolate(SpatialCoordinate(mesh)) dat = u.at(points.dat.data_ro)
I get the following error:
File "lin_dyn_unicell.py", line 117, in <module> dat = u.at(points.dat.data_ro) File "/Users/gidel/Documents/Firedrake/firedrake/src/firedrake/firedrake/function.py", line 585, in at raise PointNotInDomainError(self.function_space().mesh(), points[i].reshape(-1)) firedrake.function.PointNotInDomainError: domain <Mesh #1> does not contain point [200. 66.66666667]
How is that possible while both function spaces are defined on the same mesh?
You have run into the delightful world of computational geometry. Unfortunately points might not be determined to be in any cell if they fall close to cell boundaries. You can reduce the tightness of the tolerance used to check if a point is in a cell by saying: data = u.at(points, tolerance=1e-6) Thanks, Lawrence _______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
participants (3)
-
Floriane Gidel [RPG]
-
Ham, David A
-
Lawrence Mitchell