This is a read only copy of the old FEniCS QA forum. Please visit the new QA forum to ask questions

dofmap and nodal solutions

+4 votes

Hi,
I'm looking for a function that relates the nodal solutions (I guess should be in RAM before projection) with the dof index. The final purpose is to try to make a proper drawing 2D DG.

I show you a simplified retail of code to focus the problem.

from dolfin import *
mesh=Mesh(Rectangle(-2, -2,2,2),2)
f=conditional(lt(Expression('x[0]'),0 ),5.0, conditional(lt(Expression('x[1]'),0 ),3.0,0))

# only for more realistic visualization of f
#V = FunctionSpace(mesh, 'DG', 0) 

V = FunctionSpace(mesh, 'DG', 2)
projection=project (f,V)

Here, I relate the nodal coordinates with the dof index and I need the code to relate the projection values with dof index

dofmap = V.dofmap()
for cell in cells(mesh):
    node_coordinates = dofmap.tabulate_coordinates(cell)
    dof_index = dofmap.cell_dofs(cell.index())
    print 'node coordinates:'
    print node_coordinates
    print 'dofs:'
    print dof_index
#plot (projection)
#interactive()
#file2 = File("./sol.xml")
#file2 << projection

I do not like the only solution I've found: export to xml and read from there the 3 and 5 column.. Maybe I'm wrong but I think the data before performing the projection are somewhere ... but I can not access them ... :-(

Retail of xml code..

<dof index="0" value="3,0000000000000138" cell_index="0" cell_dof_index="0"/>
<dof index="1" value="3,0000000000000164" cell_index="0" cell_dof_index="1"/>
<dof index="2" value="3,0000000000000053" cell_index="0" cell_dof_index="2"/>

Javi

asked Jul 4, 2014 by ajaome FEniCS Novice (450 points)
edited Jul 4, 2014 by ajaome

1 Answer

+9 votes
 
Best answer

Hi, see if this helps

from dolfin import *

mesh = UnitSquareMesh(2, 2)

f = Expression('sin(pi*x[0])*x[1]')
V = FunctionSpace(mesh, 'DG', 2)

f_proj = project(f, V)

F = f_proj.vector().array()
X = V.dofmap().tabulate_all_coordinates(mesh)
X.resize((V.dim(), 2))

print 'dof index | dof coordinate |  dof value'
for i, (x, v) in enumerate(zip(X, F)):
    print i, x, v 
answered Jul 4, 2014 by MiroK FEniCS Expert (80,920 points)
selected Jul 4, 2014 by ajaome

Very clear and helpful, thank you very much. I have to study more fenics and python... but it goes slow and hard.

...