dolfin team mailing list archive

[Martin Sandve Alnæs <martinal@xxxxxxxxx>]

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 <martinal@xxxxxxxxx> 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 <745646@xxxxxxxxxxxxxxxxxx>:
>> 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@xxxxxxxxxxxxxxxxxxx
>> Unsubscribe : https://launchpad.net/~dolfin
>> More help   : https://help.launchpad.net/ListHelp
>>
>