Extract solution at a set of nodes

Question

Dear all,

After solving a problem over a mesh, I need to extract the solution at a set of nodes—basically I need to build a small table that contains coordinates and solution at the nodes. I build my mesh using GMSH, and the set of nodes belongs to a "Physical Line". I'm able to import this Physical Line as a facet mesh function, but then I can seem to figure out of to extract the solution at the nodes. Any help is appreciated. Thanks!

Martin

MiroK · Answer 1 · Dec 9, 2016

Hi, in serial consider the following

from dolfin import *
import numpy as np

mesh = UnitSquareMesh(64, 64)
facet_f = FacetFunction('size_t', mesh, 0)
CompiledSubDomain('near(x[0], x[1])').mark(facet_f, 1)

V = FunctionSpace(mesh, 'CG', 1)
v = interpolate(Expression('x[0]+x[1]', degree=1), V)

mesh.init(1, 0)
# Unique vertices of the line (by indices)
vs = list(set(sum((f.entities(0).tolist() for f in SubsetIterator(facet_f, 1)), [])))
# Get degrees of freedom associated to vertices (by indices)
v2d = vertex_to_dof_map(V)
d = v2d[vs]
# Dof values
values = v.vector().array()[d]
# The values should really be the same as evaluating v at vertices (by coordinates)
x = mesh.coordinates().reshape((-1, 2))
# Collect values at points that are on line 
values0 = np.array([v(xi) for xi in x[vs]])

print np.linalg.norm(values0 - values)

meron · Answer 2 · Dec 8, 2016

Do you mean solution.vector().array() ?

answered Dec 8, 2016 by meron FEniCS User (2,340 points)

Thanks for your help! However, this gives me the solution for dofs. How can I extract only the solution for the nodes that belong to a given line? Thanks again. Martin

commented Dec 9, 2016 by Martin Genet FEniCS User (1,460 points)

Maybe have a look at this: https://answers.launchpad.net/dolfin/+question/224725
Else, I wouldn't know either, sorry.
Cheers, Meron

commented Dec 9, 2016 by meron FEniCS User (2,340 points)

hirshikesh · Answer 3 · Dec 9, 2016

        coor = mesh.coordinates()            
        for i in range(len(u_array)):
           if coor[i][1]==0.0:   # fixed line y = 0.0
              print 'u(%8g,%8g) = %g' % (coor[i][0], coor[i][1], u_array[i])

answered Dec 9, 2016 by hirshikesh FEniCS User (1,890 points)

Thanks for your help! However, my line is actually a complex path. Any idea how to extract the solution at its nodes? Plus, I am not sure index i of mesh.coordinates() and solution.vector().array() correspond, see https://fenicsproject.org/qa/2715/coordinates-u_nodal_values-using-numerical-source-function. Anyway, thanks again. Martin

commented Dec 9, 2016 by Martin Genet FEniCS User (1,460 points)

Extract solution at a set of nodes

3 Answers