dolfin team mailing list archive

[Garth Wells <785874@xxxxxxxxxxxxxxxxxx>]

The NaN is not caused by FFC autodetecting the degree.

The numerical error in some cases is due to FFC not increasing the
degree when the spatial coordinate is used as a coefficient function.

-- 
You received this bug notification because you are a member of DOLFIN
Team, which is subscribed to DOLFIN.
https://bugs.launchpad.net/bugs/785874

Title:
  Projection of x is not accurate

Status in DOLFIN:
  New

Bug description:
  I've tested that projecting x works without the scaling bug that was
  just fixed, using dimensions 1,2,3 and both DG and CG from 0 to 3
  degrees. I print the max and min values of the vector of the
  projection function, and the values are _close_ to 0 and 1 but not to
  machine precision. The script is below.

  There's up to 2.7% error in the 3D case. Is the projection form
  integrated accurately enough? All but the DG0 function space should be
  capable of representing x exactly. Not sure if this is a dolfin or ffc
  bug.

  
  from dolfin import *

  def mcx(dim):
      if dim == 1:
          mesh = UnitInterval(20)
          cell = interval
          x = cell.x
      if dim == 2:
          mesh = UnitSquare(20, 20)
          cell = triangle
          x = cell.x[0]
      if dim == 3:
          mesh = UnitCube(20, 20, 20)
          cell = tetrahedron
          x = cell.x[0]
      return mesh, cell, x

  for dim in range(1, 4):
      mesh, cell, x = mcx(dim)
      minval, maxval = 1.0, 0.0
      #print dim, "DG"
      for degree in range(3):
          V = FunctionSpace(mesh, "DG", degree)
          u = project(x, V)
          #print dim, degree, u.vector().min(), u.vector().max()
          minval = min(minval, u.vector().min())
          maxval = max(maxval, u.vector().max())
      #print dim, "CG"
      for degree in range(1, 4):
          V = FunctionSpace(mesh, "CG", degree)
          u = project(x, V)
          #print dim, degree, u.vector().min(), u.vector().max()
          minval = min(minval, u.vector().min())
          maxval = max(maxval, u.vector().max())
      print minval, maxval