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
