← Back to team overview

ufl team mailing list archive

Re: [FFC-dev] bug with tensor fields

 

On Tuesday 22 September 2009 10:03:12 kent-and@xxxxxxxxx wrote:
> > On Tuesday 22 September 2009 09:49:38 kent-and@xxxxxxxxx wrote:
> >> The following variational forms should produce the same matrices,
> >>
> >> Ai = Function(T, ( ("1.0", "0.0"), ("0.0", "1.0")) )
> >> a1 = dot(dot(Ai, curl(v)), curl(u))*dx
> >> a2 = dot(curl(v), curl(u))*dx
> >>
> >> But it doesn't.
> >>
> >> I am not sure whether the bug is in Dolfin, FFC or UFL.
> >>
> >> The following test script can be used to test the problem.
> >>
> >> from dolfin import *
> >>
> >> N = 2
> >> mesh = UnitCube(N,N,N)
> >> V = FunctionSpace(mesh, "N1curl", 1)
> >> T = TensorFunctionSpace(mesh, "CG", 1)
> >>
> >> u = TrialFunction(V)
> >> v = TestFunction(V)
> >>
> >> Ai = Function(T, ( ("1.0", "0.0"), ("0.0", "1.0")) )
> >
> > This creates a 2x2 tensor function, do you not want a 3x3 tensor?
> >
> > I have plans on fixing rigorous tests for value shapes for function
> > expressions. Which are automatically checked with the value shape of the
> > FunctionSpace. This could prevent your error and some nasty segfaults...
> >
> > Johan
> 
> That would be great. I think a simple check on the value shape and
> the tuples is sufficient.

Yes something like that. What makes it not straightforwardly is that a user 
can define his own Function class:

  class MyFunc(Function):
     cpparg = (("2","3"),("3","4"))

Then I need to create some value_shape function that is used when this 
function is instantiated with a FunctionSpace.

It is not nessesary to have this for pure Python expressions, like:

  class MyFunc(Function):
       def eval(values,x):
          .....

Here will the callback typemap produce the correct value shape on the NumPy 
array that is passed to the user, with all its index checks.

Johan

> Kent
> 
> >> a1 = dot(dot(Ai, curl(v)), curl(u))*dx
> >> a2 = dot(curl(v), curl(u))*dx
> >>
> >> A1 = assemble(a1)
> >> A2 = assemble(a2)
> >>
> >>
> >> file = File("A1.m")
> >> file <<A1
> >> file = File("A2.m")
> >> file <<A2
> >>
> >> file_str = """
> >> A1;
> >> A1 = A;
> >> clear A;
> >> A2;
> >> A2 = A;
> >> clear A;
> >>
> >> max(max(A1-A2))
> >> """
> >>
> >> f = open("test_tensor.m", 'w')
> >> f.write(file_str)
> >> f.close()
> >>
> >> import os
> >> os.system("octave test_tensor.m")
> >>
> >>
> >>
> >> _______________________________________________
> >> FFC-dev mailing list
> >> FFC-dev@xxxxxxxxxx
> >> http://www.fenics.org/mailman/listinfo/ffc-dev
> 


References