Dear Koki and Lawrence, I had reached out to you On Feb 23rd regarding an issue I am facing with regards to extracting results at certain points. I have also started issue #2359 in the github page. I apologize to bother you like this but I would be very grateful if you could kindly advise on how I may resolve the issue. -- Abhishek Venketeswaran NETL Research Associate - ORISE National Energy Technology Laboratory Department of Energy abhishek.venketeswaran@netl.doe.gov Work: 412-386-4833 Mobile: 716-507-7890 www.netl.doe.gov<http://www.netl.doe.gov> [cid:311ae572-2f97-43a1-a610-3dbcd81a2782] ________________________________ From: Venketeswaran, Abhishek (FN) Sent: Wednesday, February 23, 2022 3:22 PM To: Sagiyama, Koki <k.sagiyama@imperial.ac.uk>; Lawrence Mitchell <wence@gmx.li>; firedrake@imperial.ac.uk <firedrake@imperial.ac.uk> Subject: Error with Point Evaluation and VertexOnlyMesh Dear Koki and Lawrence, I am having trouble extracting my results at certain points on my mesh. I understand that there are 2 methods to do this (i) at() function and (ii) VertexOnlyMesh. Both dont work for me. While the former raises Point Not in Domain Error, the latter raises TypeError: an integer is required. Here is a minimal example of trying both methods out: from firedrake import * from firedrake.petsc import PETSc from firedrake.mesh import VertexOnlyMesh from collections import OrderedDict import numpy as np mesh = UnitSquareMesh(10, 10) V = FunctionSpace(mesh, "CG", 1) v = TestFunction(V) u = TrialFunction(V) a = inner(grad(u), grad(v)) * dx L = Constant(0) * v * dx bc1 = DirichletBC(V, 0.0, 1) bc2 = DirichletBC(V, 1.0, 2) sol = Function(V) problem1 = LinearVariationalProblem(a, L, sol, bcs=[bc1, bc2], constant_jacobian=True) solver1 = LinearVariationalSolver(problem1, solver_parameters={"ksp_type": "cg", "pc_type": "jacobi"}) solver1.solve() print(sol.dat.data_ro) crd=[0.5256,0.5051] num_points = 25 δs = np.linspace(0.0, 1.0, num_points + 1) X, Y = np.meshgrid(δs, δs) xs = np.vstack((X.flatten(), Y.flatten())).T vm=1 pe=1 #Test VertexOnlyMesh Generation if vm: point_cloud = VertexOnlyMesh(mesh, xs) #Test Point Evaluation of Function if pe: print(sol.at(crd)) Kindly let me know where I should start looking to fix this problem. -- Abhishek Venketeswaran NETL Research Associate - ORISE National Energy Technology Laboratory Department of Energy abhishek.venketeswaran@netl.doe.gov Work: 412-386-4833 Mobile: 716-507-7890 www.netl.doe.gov<http://www.netl.doe.gov> [cid:b3222045-e3c2-4f40-9a6f-b6c662401b72]
