Visualisation of dofs when projecting function

[Jørgen Dokken]

When projecting a simple expression, in this case x*y onto a CG space the visualization of dofs disappears when the projected value is a small negative number,

f(x,y)<0, abs(f(x,y))<10e-17.

I'm using DOLFIN 2016.2.0.dev0. This bug is not in DOLFIN 1.6

from dolfin import *
mesh = UnitSquareMesh(15,15)
V = FunctionSpace(mesh,'CG',1)
f = Expression('x[0]*x[1]')
q = project(f,V)

[Jan Blechta]

Following snippet tries to explain the behaviour

from __future__ import print_function
from dolfin import *

mesh = UnitSquareMesh(9, 3)
V = FunctionSpace(mesh, 'CG', 1)

print("Projection from P1 to P1 would ideally be identity mapping")
f = Expression('x[0]*x[1]', degree=1)

print("Negative value is numerical roundoff in asembly, LU solve")
q = project(f, V, solver_type='lu')
print(q.vector().min(), q.vector().max())

print("Negative value is iterative error in CG approximation")
q = project(f, V, solver_type='cg')
print(q.vector().min(), q.vector().max())
print()

print("Projection from P2 to P1 is not subject to maximum principle,")
print("negative value expected")
f = Expression('x[0]*x[1]', degree=2)
q = project(f, V)
print(q.vector().min(), q.vector().max())
print()

print("Lagrange interpolation of non-negative function is non-negative")
f = Expression('x[0]*x[1]', degree=2)
q = interpolate(f, V)
print(q.vector().min(), q.vector().max())

Output is

Projection from P1 to P1 would ideally be identity mapping
Negative value is numerical roundoff in asembly, LU solve
-4.27028924952e-17 1.0
Negative value is iterative error in CG approximation
-3.13058291719e-07 1.00000007143

Projection from P2 to P1 is not subject to maximum principle,
negative value expected
-0.00451511219686 0.995484887803

Lagrange interpolation of non-negative function is non-negative
0.0 1.0

[Jørgen Dokken]

But why does the visualization remove these nodes?

[Jan Blechta]

Could you post complete code to reproduce the problem?

[Jørgen Dokken]

The only thing added to the code above is

plot(q)

[Jan Blechta]

Can't reproduce. Is the vector correct? What's the output of print(q.vector().array())?

[Jørgen Dokken]

from dolfin import *
import numpy as np

mesh = UnitSquareMesh(10,10)
V = FunctionSpace(mesh,'CG',1)
f = Expression('x[0]*x[1]')
q = project(f,V)
plot(q,title='Plot of x*y')
interactive()
t= q.vector().array()
print t

Gives the following output:

[ -4.34353546e-18   8.68707093e-18   1.00000000e-01   4.07119023e-18
   9.00000000e-02   2.00000000e-01   7.83378043e-18   8.00000000e-02
   1.80000000e-01   3.00000000e-01   7.02331485e-18   7.00000000e-02
   1.60000000e-01   2.70000000e-01   4.00000000e-01  -9.61382395e-18
   6.00000000e-02   1.40000000e-01   2.40000000e-01   3.60000000e-01
   5.00000000e-01  -1.35525272e-18   5.00000000e-02   1.20000000e-01
   2.10000000e-01   3.20000000e-01   4.50000000e-01   6.00000000e-01
   8.13151629e-18   4.00000000e-02   1.00000000e-01   1.80000000e-01
   2.80000000e-01   4.00000000e-01   5.40000000e-01   7.00000000e-01
   8.52951478e-18   3.00000000e-02   8.00000000e-02   1.50000000e-01
   2.40000000e-01   3.50000000e-01   4.80000000e-01   6.30000000e-01
   8.00000000e-01  -2.38799094e-18   2.00000000e-02   6.00000000e-02
   1.20000000e-01   2.00000000e-01   3.00000000e-01   4.20000000e-01
   5.60000000e-01   7.20000000e-01   9.00000000e-01   3.04931861e-18
   1.00000000e-02   4.00000000e-02   9.00000000e-02   1.60000000e-01
   2.50000000e-01   3.60000000e-01   4.90000000e-01   6.40000000e-01
   8.10000000e-01   1.00000000e+00  -5.42101086e-18   2.00000000e-02
   6.00000000e-02   1.20000000e-01   2.00000000e-01   3.00000000e-01
   4.20000000e-01   5.60000000e-01   7.20000000e-01   9.00000000e-01
  -8.44018659e-19   3.00000000e-02   8.00000000e-02   1.50000000e-01
   2.40000000e-01   3.50000000e-01   4.80000000e-01   6.30000000e-01
   8.00000000e-01  -6.63157466e-18   4.00000000e-02   1.00000000e-01
   1.80000000e-01   2.80000000e-01   4.00000000e-01   5.40000000e-01
   7.00000000e-01  -1.08420217e-17   5.00000000e-02   1.20000000e-01
   2.10000000e-01   3.20000000e-01   4.50000000e-01   6.00000000e-01
   1.63053842e-17   6.00000000e-02   1.40000000e-01   2.40000000e-01
   3.60000000e-01   5.00000000e-01   2.71050543e-18   7.00000000e-02
   1.60000000e-01   2.70000000e-01   4.00000000e-01  -4.13352078e-17
   8.00000000e-02   1.80000000e-01   3.00000000e-01  -3.71964745e-18
   9.00000000e-02   2.00000000e-01   1.02174140e-17   1.00000000e-01
  -8.13151629e-18]

and following plot:

[Jan Blechta]

Can't reproduce the problem.

[Jørgen Dokken]

Really strange, When I ran this I was using the latest pull of dolfin: 2016.2.0.dev0. Do you have any idea of it is the version of some dependencies that may cause this issue?

[Jan Blechta]

Looks like a VTK issue. Does your build has a general problem plotting negative values?

[Jørgen Dokken]

I've never experienced issues with negative numbers earlier. This is my current VTK-setup

dokken@dokken ~/Desktop $ dpkg -l | grep vtk
ii  libvtk5.10                                 5.10.1+dfsg-2.1build1                                       amd64        Visualization Toolkit - A high level 3D visualization library - runtime
ii  libvtk6-java                               6.2.0+dfsg1-10build1                                        amd64        Visualization Toolkit - A high level 3D visualization library - java
ii  libvtk6.2                                  6.2.0+dfsg1-10build1                                        amd64        VTK libraries
ii  libvtk6.2-qt                               6.2.0+dfsg1-10build1                                        amd64        VTK libraries, Qt files
ii  python-vtk                                 5.10.1+dfsg-2.1build1                                       amd64        Python bindings for VTK
ii  tcl-vtk                                    5.10.1+dfsg-2.1build1                                       amd64        Tcl bindings for VTK
ii  tcl-vtk6                                   6.2.0+dfsg1-10build1                                        amd64        Tcl bindings for VTK
ii  vtk6                                       6.2.0+dfsg1-10build1                                        amd64        Binaries for VTK6