On 8 June 2011 12:11, Martin Sandve Alnæs <marti...@simula.no> wrote:
> Done and checked in. If someone updates FFC to support this, we can
> hopefully close this bug.
I'm not sure this is enough to handle the bug. If you look at the
output of printing M in the example code I posted you'll see that the
list tensor contains component '7'. This is what you'll get from
calling self.component(), but the
TT.symmetry() only contains:
{(2,): (1,), (6,): (5,)}
Is there some function we need to call first to map the component '7'
--> '6', before looking at symmetries to map '6' --> '5'?
Doing so will get us into trouble with mapping '3' --> '2' since
symmetry will map that to '1'.
The TT element has 2 x 3 sub elements due to symmetry.
Kristian
> Martin
>
>
> On 8 June 2011 11:17, Martin Sandve Alnæs <marti...@simula.no> wrote:
>> Update: Of course the expand_indices algorithm does not handle your
>> example where the tensor element is part of a mixed element. I'll
>> consider adding symmetry as a property of all elements to fix this,
>> that is much less work than adding it to all expressions.
>>
>> Martin
>>
>> On 8 June 2011 10:58, Martin Sandve Alnæs <marti...@simula.no> wrote:
>>> I checked out what UFL does. If you look at expand_indices.py,
>>> you see that this algorithm does handle symmetries. However,
>>> that algorithm is not (and should not be) used by all form compiler
>>> variants. Here's a simplified version of the code:
>>>
>>> def form_argument(self, x):
>>> if x.shape() == ():
>>> return x
>>> else:
>>> e = x.element()
>>>
>>> # Get component
>>> c = self.component()
>>>
>>> # Map it through the symmetry mapping
>>> if isinstance(e, TensorElement):
>>> s = e.symmetry() or {}
>>> c = s.get(c, c)
>>>
>>> return x[c]
>>>
>>>
>>> I want to make two points from this code:
>>>
>>> 1) The symmetry mapping is currently only available directly from a
>>> TensorElement. Extending symmetries to be a property of all ufl
>>> expressions could be quite a bit of work. If this is a wanted feature,
>>> make a blueprint.
>>>
>>> 2) When you do have the symmetry mapping at hand, the mapping is
>>> trivial. It should be a quick fix in FFC, just apply the symmetry
>>> mapping like above everywhere you get the component of a form
>>> argument.
>>>
>>> Martin
>>>
>>> 2011/4/5 Kristian B. Ølgaard <745...@bugs.launchpad.net>:
>>>> The following simple form also fails:
>>>>
>>>> T1 = TensorElement('CG', triangle, 1, symmetry=True)
>>>> T2 = TensorElement('CG', triangle, 1, symmetry=True)
>>>> TT = T1*T2
>>>> P, Q = Coefficients(TT)
>>>> M = inner(P, Q)*dx
>>>> print M
>>>>
>>>> Printing 'M' results in:
>>>>
>>>> { ([[
>>>> [ (w_0)[0], (w_0)[1] ],
>>>> [ (w_0)[1], (w_0)[3] ]
>>>> ]]) : ([[
>>>> [ (w_0)[4], (w_0)[5] ],
>>>> [ (w_0)[5], (w_0)[7] ]
>>>> ]]) } * dx0
>>>>
>>>> which illustrates the problem in FFC for both tensor and quadrature
>>>> representations.
>>>> The component '7' in the ListTensor does not exist in the MixedElement of
>>>> FFC which, due to symmetry, only contain 6 'unique' subelements
>>>> (components).
>>>> One could argue that UFL should keep track of symmetry when creating the
>>>> indices of list tensors such that it maps 3->2, 4->3, 5->4 and 7->6.
>>>>
>>>> --
>>>> You received this bug notification because you are a member of DOLFIN
>>>> Team, which is subscribed to DOLFIN.
>>>> https://bugs.launchpad.net/bugs/745646
>>>>
>>>> Title:
>>>> Problem with assemble() with MixedFunctionSpace of symmetric
>>>> TensorFunctionSpaces
>>>>
>>>> Status in DOLFIN:
>>>> New
>>>>
>>>> Bug description:
>>>> from dolfin import *
>>>>
>>>> mesh = Rectangle(0,0,1,1,10,10)
>>>> dim = 2 #assume 2D
>>>> symm = dict(((i,j), (j,i))
>>>> for i in range(dim) for j in range(dim) if i > j )
>>>>
>>>> T1 = TensorFunctionSpace(mesh, 'CG', 1, symmetry=symm)
>>>> T2 = TensorFunctionSpace(mesh, 'CG', 1, symmetry=symm)
>>>>
>>>> TT = T1*T2
>>>>
>>>> R, S = TrialFunctions(TT)
>>>> v_R, v_S = TestFunctions(TT)
>>>>
>>>> P = Function(T1)
>>>> Q = Function(T2)
>>>>
>>>> Fr = inner(R,v_R)*dx+ inner(P,v_R)*dx
>>>>
>>>> Fs = inner( S, v_S )*dx + inner( Q, v_S )*dx
>>>>
>>>> F = Fr + Fs
>>>>
>>>> A = lhs(F)
>>>> B = rhs(F)
>>>>
>>>> a = Matrix()
>>>> a = assemble(A, tensor=a)
>>>> b = Vector()
>>>> b = assemble(B, tensor=b)
>>>>
>>>> I get the following error:
>>>>
>>>> Calling FFC just-in-time (JIT) compiler, this may take some time.
>>>> Traceback (most recent call last):
>>>> File "tensortest2.py", line 28, in <module>
>>>> a = assemble(A, tensor=a)
>>>> File "/usr/lib/python2.6/dist-packages/dolfin/fem/assemble.py", line
>>>> 100, in assemble
>>>> common_cell=common_cell)
>>>> File "/usr/lib/python2.6/dist-packages/dolfin/fem/form.py", line 34, in
>>>> __init__
>>>> (self._compiled_form, module, self.form_data) = jit(form,
>>>> form_compiler_parameters, common_cell)
>>>> File "/usr/lib/python2.6/dist-packages/dolfin/compilemodules/jit.py",
>>>> line 47, in mpi_jit
>>>> return local_jit(*args, **kwargs)
>>>> File "/usr/lib/python2.6/dist-packages/dolfin/compilemodules/jit.py",
>>>> line 114, in jit
>>>> return jit_compile(form, parameters=p, common_cell=common_cell)
>>>> File "/usr/lib/python2.6/dist-packages/ffc/jitcompiler.py", line 64, in
>>>> jit
>>>> return jit_form(object, parameters, common_cell)
>>>> File "/usr/lib/python2.6/dist-packages/ffc/jitcompiler.py", line 122,
>>>> in jit_form
>>>> compile_form(preprocessed_form, prefix=jit_object.signature(),
>>>> parameters=parameters)
>>>> File "/usr/lib/python2.6/dist-packages/ffc/compiler.py", line 140, in
>>>> compile_form
>>>> ir = compute_ir(analysis, parameters)
>>>> File "/usr/lib/python2.6/dist-packages/ffc/representation.py", line 66,
>>>> in compute_ir
>>>> irs = [_compute_integral_ir(f, i, parameters) for (i, f) in
>>>> enumerate(forms)]
>>>> File "/usr/lib/python2.6/dist-packages/ffc/representation.py", line
>>>> 186, in _compute_integral_ir
>>>> parameters)
>>>> File
>>>> "/usr/lib/python2.6/dist-packages/ffc/tensor/tensorrepresentation.py",
>>>> line 59, in compute_integral_ir
>>>> ir["AK"] = _compute_terms(monomial_form, None, None, domain_type,
>>>> quadrature_degree)
>>>> File
>>>> "/usr/lib/python2.6/dist-packages/ffc/tensor/tensorrepresentation.py",
>>>> line 98, in _compute_terms
>>>> quadrature_degree)
>>>> File "/usr/lib/python2.6/dist-packages/ffc/tensor/referencetensor.py",
>>>> line 28, in __init__
>>>> self.A0 = integrate(monomial, domain_type, facet0, facet1,
>>>> quadrature_order)
>>>> File
>>>> "/usr/lib/python2.6/dist-packages/ffc/tensor/monomialintegration.py", line
>>>> 50, in integrate
>>>> psis = [_compute_psi(v, table, len(points), domain_type) for v in
>>>> monomial.arguments]
>>>> File
>>>> "/usr/lib/python2.6/dist-packages/ffc/tensor/monomialintegration.py", line
>>>> 169, in _compute_psi
>>>> Psi[component][tuple(dlist)] = etable[dtuple][:,
>>>> cindex[0].index_range[component], :]
>>>>
>>>> _______________________________________________
>>>> Mailing list: https://launchpad.net/~dolfin
>>>> Post to : dolfin@lists.launchpad.net
>>>> Unsubscribe : https://launchpad.net/~dolfin
>>>> More help : https://help.launchpad.net/ListHelp
>>>>
>>>
>>
>
> _______________________________________________
> Mailing list: https://launchpad.net/~dolfin
> Post to : dolfin@lists.launchpad.net
> Unsubscribe : https://launchpad.net/~dolfin
> More help : https://help.launchpad.net/ListHelp
>
_______________________________________________
Mailing list: https://launchpad.net/~dolfin
Post to : dolfin@lists.launchpad.net
Unsubscribe : https://launchpad.net/~dolfin
More help : https://help.launchpad.net/ListHelp
