Spatial coordinates in form are scaled by 2/3 in 2D

[Martin Sandve Alnæs]

Here's a program that plots y=x, y=2x, y=3x in 1D and 2D, giving correct plots in 1D but scaled by 2/3 in 2D:

from dolfin import *

n = 30
factor = 1.0

#cell = interval
cell = triangle

d = cell.d
if d == 1:
    mesh = UnitInterval(n)
elif d == 2:
    mesh = UnitSquare(n, n)
    #factor = 3.0/2.0 # Enable this to adjust for bug and get correct plots

# Asserting that x indeed goes from 0 to 1
assert mesh.coordinates()[:].max() == 1.0
assert mesh.coordinates()[:].min() == 0.0

# Hack for bug in ffc, fix proposed for merging into dev branch
x = cell.x if d == 1 else cell.x[0]

V = FunctionSpace(mesh, "CG", 1)

for k in range(1, 4):
    u = factor*k*x
    print str(u)
    uh = project(u, V)
    plot(uh, vmin=-1.0, vmax=1.0)

interactive()

[Garth Wells (garth-wells)]

This test is a bit cryptic. Can you make it easy for us and explain what the expected result is (and how to check it), and what the actual result is (and how to check it).

[Martin Sandve Alnæs (martinal)]

Written on the bus and untested:

mesh = UnitSquare(30,30)
V = FunctionSpace(mesh, "CG", 1)
y = project(triangle.x[0], V)
print y.vector().max() # should be 1.0, is 0.67

Martin
Den 19. mai 2011 18.10 skrev "Garth Wells" <email address hidden>
følgende:
> This test is a bit cryptic. Can you make it easy for us and explain what
> the expected result is (and how to check it), and what the actual result
> is (and how to check it).
>
> --
> You received this bug notification because you are a direct subscriber
> of the bug.
> https://bugs.launchpad.net/bugs/785244
>
> Title:
> Spatial coordinates in form are scaled by 2/3 in 2D
>
> Status in DOLFIN:
> New
>
> Bug description:
> Here's a program that plots y=x, y=2x, y=3x in 1D and 2D, giving
> correct plots in 1D but scaled by 2/3 in 2D:
>
>
> from dolfin import *
>
> n = 30
> factor = 1.0
>
> #cell = interval
> cell = triangle
>
> d = cell.d
> if d == 1:
> mesh = UnitInterval(n)
> elif d == 2:
> mesh = UnitSquare(n, n)
> #factor = 3.0/2.0 # Enable this to adjust for bug and get correct plots
>
> # Asserting that x indeed goes from 0 to 1
> assert mesh.coordinates()[:].max() == 1.0
> assert mesh.coordinates()[:].min() == 0.0
>
> # Hack for bug in ffc, fix proposed for merging into dev branch
> x = cell.x if d == 1 else cell.x[0]
>
> V = FunctionSpace(mesh, "CG", 1)
>
> for k in range(1, 4):
> u = factor*k*x
> print str(u)
> uh = project(u, V)
> plot(uh, vmin=-1.0, vmax=1.0)
>
> interactive()
>
> To unsubscribe from this bug, go to:
> https://bugs.launchpad.net/dolfin/+bug/785244/+subscribe

[Martin Sandve Alnæs (martinal)]

New and improved with comments explaining:

from dolfin import *
mesh = UnitSquare(30,30)
V = FunctionSpace(mesh, "CG", 1)
x = triangle.x

# Compute f(x,y) = x on Omega = [0,1]^2
f = project(x[0], V)

# Assert that max f = max x on Omega = 1
fmax = f.vector().max()
print fmax # Prints 0.67, not 1.0
assert abs(fmax - 1.0) < 1e-8 # should pass, fails

On 19 May 2011 18:21, Martin Sandve Alnæs <email address hidden> wrote:
> Written on the bus and untested:
>
> mesh = UnitSquare(30,30)
>
> V = FunctionSpace(mesh, "CG", 1)
> y = project(triangle.x[0], V)
> print y.vector().max() # should be 1.0, is 0.67
>
> Martin
>
> Den 19. mai 2011 18.10 skrev "Garth Wells" <email address hidden>
> følgende:
>> This test is a bit cryptic. Can you make it easy for us and explain what
>> the expected result is (and how to check it), and what the actual result
>> is (and how to check it).
>>
>> --
>> You received this bug notification because you are a direct subscriber
>> of the bug.
>> https://bugs.launchpad.net/bugs/785244
>>
>> Title:
>> Spatial coordinates in form are scaled by 2/3 in 2D
>>
>> Status in DOLFIN:
>> New
>>
>> Bug description:
>> Here's a program that plots y=x, y=2x, y=3x in 1D and 2D, giving
>> correct plots in 1D but scaled by 2/3 in 2D:
>>
>>
>> from dolfin import *
>>
>> n = 30
>> factor = 1.0
>>
>> #cell = interval
>> cell = triangle
>>
>> d = cell.d
>> if d == 1:
>> mesh = UnitInterval(n)
>> elif d == 2:
>> mesh = UnitSquare(n, n)
>> #factor = 3.0/2.0 # Enable this to adjust for bug and get correct plots
>>
>> # Asserting that x indeed goes from 0 to 1
>> assert mesh.coordinates()[:].max() == 1.0
>> assert mesh.coordinates()[:].min() == 0.0
>>
>> # Hack for bug in ffc, fix proposed for merging into dev branch
>> x = cell.x if d == 1 else cell.x[0]
>>
>> V = FunctionSpace(mesh, "CG", 1)
>>
>> for k in range(1, 4):
>> u = factor*k*x
>> print str(u)
>> uh = project(u, V)
>> plot(uh, vmin=-1.0, vmax=1.0)
>>
>> interactive()
>>
>> To unsubscribe from this bug, go to:
>> https://bugs.launchpad.net/dolfin/+bug/785244/+subscribe
>

[Garth Wells (garth-wells)]

Below is a test program that digs a bit deeper. Looks like an FFC issue. Only seems to impact zero and first order elements.

from dolfin import *
mesh = UnitSquare(30, 30)
V = FunctionSpace(mesh, "DG", 0)
v = TestFunction(V)

x = triangle.x
vol = Cell(mesh, 0).volume()
L = x[0]*v*dx
b = assemble(L)
fmax = b.max()/Cell(mesh, 0).volume()

# This should be close to 1.0, but it's close to 2.0
print fmax

[Garth Wells (garth-wells)]

Bug is an index error in the array

  FEA1_f0

in the code generated by FFC. Will try to fix.

[Garth Wells (garth-wells)]

Fixed in FFC (http://bazaar.launchpad.net/~ffc-core/ffc/main/revision/1640)