Subject: Re: [FFC-dev] [DOLFIN-dev] Mixed/Quadrature element

From: Kristian Oelgaard <[email-masked]>
Date: Wed, 07 May 2008 11:32:21 +0200

Quoting Jake Ostien <[email-masked]>:

> Hi,
> 
> I am trying to define a mixed element with one part linear Lagrange 
> elements and one part QuadratureElement.  The error message that I get 
> tells me that I can't take derivatives of a QuadratureElement, but I 
> don't think that I am.  Please see the attached form.  Is this something 
> that I should be able to do?

I think that you should be able to do this. Somehow the compiler applies the
derivative to the mixed element (which is made up of two mixed elements) where
it should only take the derivative of the first of the mixed elements.

Printing the basisfunctions in monomialintegration.py __init_table() I get:

[((d/dXa0[0, 1, 2])vi0[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16,
17, 18, 19, 20][b0[0, 1, 2]]), ((d/dXa1[0, 1, 2])vi1[0, 1, 2, 3, 4, 5, 6, 7, 8,
9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20][b0[0, 1, 2]])]

Then, when initialising the table, the compiler will ask the finite elements to
tabulate derivatives up to order 1, but since derivatives are not defined for
QuadratureElements it breaks. I'm not sure where the best place to fix this
would be...

Kristian

> Jake

elu  = VectorElement("Lagrange", "tetrahedron", 1)
elHp = VectorQuadratureElement("tetrahedron", 1, 9)

sgp = MixedElement([elu, elHp])

(w, V) = TestFunctions(sgp)
(u, Hp) = TrialFunctions(sgp)

def Mat(A):
    return [ [A[0], A[1], A[2]],    \
             [A[3], A[4], A[5]],    \
             [A[6], A[7], A[8]] ]

def sym(v):
    return 0.5*(v + transp(v))

mu = 0.3
lmbda = 10.0

def sigma(A):
    return mult(2.0*mu,A) + mult(lmbda,mult(trace(A),Identity(3)))

a = dot(grad(w), sigma(sym(grad(u))))*dx

---

From: Anders Logg <[email-masked]>
Date: Wed, 7 May 2008 19:20:54 +0200

> QuadratureElements it breaks. I'm not sure where the best place to fix this

I suggest manually making sure that you only take derivatives of the
part which is not a quadrature element (by reformulating the form).

This (as everything else) will be fixed in UFL. :-)

-- 
Anders