# How to apply quadratic element type

Hello,

I recently found a code for the cantilever beam with uniform distributed load problem. It works well with CG element type. I want it run with quadratic element type. So I tried to edit some codes inside.
The code is here:

from future import print_function
from dolfin import *
from fenics import *
import matplotlib.pyplot as plt
import math

L = 25.
H = 1.
Nx = 250
Ny = 10
mesh = RectangleMesh(Point(0., 0.), Point(L, H), Nx, Ny, “crossed”)

def eps(v):

E = Constant(1e5)
nu = Constant(0.3)

lmbda = Enu/(1+nu)/(1-2nu)
mu = E/2/(1+nu)
lmbda = 2mulmbda/(lmbda+2*mu)

def sigma(v):
return lmbdatr(eps(v))Identity(2) + 2.0mueps(v)

rho_g = 1e-3
f = Constant((0, -rho_g))

V = VectorFunctionSpace(mesh, fe, degree = 1)
du = TrialFunction(V)
u_ = TestFunction(V)
a = inner(sigma(du), eps(u_))*dx
l = inner(f, u_)*dx

def left(x, on_boundary):
return near(x[0],0.)

bc = DirichletBC(V, Constant((0.,0.)), left)

u = Function(V, name=“Displacement”)
solve(a == l, u, bc)

plot(1e3*u, mode=“displacement”)
plt.title(“2D elastic bending model”)
plt.xlabel(“length”)
plt.ylabel(“bending”)
plt.show()

print(“Maximal deflection:”, -u(L,H/2.)[1])
print(“Beam theory deflection:”, float(3rho_gL4/2/E/H3))

Vsig = TensorFunctionSpace(mesh, “DG”, degree=0)
sig = Function(Vsig, name=“Stress”)
sig.assign(project(sigma(u), Vsig))
print(“Stress at (0,H):”, sig(0, H))

error = abs(-u(L,H/2.)[1] - float(3rho_gL4/2/E/H3))
print(error)

But the an err occurred:

runcell(0, ‘/Users/pigdtt12345/Desktop/Bill/FEniCS/Training/2D_elasticity_simple.py’)
Calling FFC just-in-time (JIT) compiler, this may take some time.
Calling FFC just-in-time (JIT) compiler, this may take some time.
Calling FFC just-in-time (JIT) compiler, this may take some time.
Symmetric part of tensor with rank != 2 is undefined.
Traceback (most recent call last):

File “/Users/pigdtt12345/Desktop/Bill/FEniCS/Training/2D_elasticity_simple.py”, line 34, in
a = inner(sigma(du), eps(u_))*dx

File “/Users/pigdtt12345/Desktop/Bill/FEniCS/Training/2D_elasticity_simple.py”, line 25, in sigma
return lmbdatr(eps(v))Identity(2) + 2.0mueps(v)

File “/Users/pigdtt12345/Desktop/Bill/FEniCS/Training/2D_elasticity_simple.py”, line 14, in eps

File “/Users/pigdtt12345/opt/anaconda3/envs/fenicsproject/lib/python3.8/site-packages/ufl/operators.py”, line 318, in sym
return Sym(A)

File “/Users/pigdtt12345/opt/anaconda3/envs/fenicsproject/lib/python3.8/site-packages/ufl/tensoralgebra.py”, line 457, in new
error(“Symmetric part of tensor with rank != 2 is undefined.”)

File “/Users/pigdtt12345/opt/anaconda3/envs/fenicsproject/lib/python3.8/site-packages/ufl/log.py”, line 172, in error
raise self._exception_type(self._format_raw(*message))

UFLException: Symmetric part of tensor with rank != 2 is undefined.

Could some please help me with it? I just want it work with “Quadrature” element type. I am sorry if I make some format error here. I am a new member with this qa forum.

If what you want is to use quadratic shape functions, I think you just need to set degree=2.

Hi,
Thanks for your help. I have tried it before.
I directly tried

V = VectorFunctionSpace(mesh, “Quadrature”, degree = 2)

but the error comes out too. The error said: Missing quad_scheme in quadrature element. So I checked online to solve how to fix it, and I found this page that has same issue.

But when I try it, it has the exception as the initial problem as “Symmetric part of tensor with rank != 2 is undefined.”

So here I am totally lost what I can do next to fix this bug. Sorry that I am beginner and this is a bit out of my knowledge.

Additionally for using quadrature elements, you might have to change the representation for the form compiler. The default is `uflacs`. This can be done using something like:

``````parameters["form_compiler"]["representation"] = 'quadrature'
# And to suppress `depecration` warnings
import warnings
``````

see here

Hello,
you are mixing “quadratic” interpolation (implictly Lagrange P2 elements) with “quadrature” elements which correspond to fields evaluated at quadrature points only and from which one cannot compute gradients and so on.
For quadratic elements, just use ` V = VectorFunctionSpace(mesh, "CG", 2)`

3 Likes

Ah! I see! Thank you so much!

Hello,

I put it at the top and the error occurs again. It said, “Symmetric part of tensor with rank != 2 is undefined.”. But I think I mixed the quadratic interpolation with quadrature elements. But thank you for your help very much too! Sorry for taking up your time.

Yeah. I just figured that you wanted to use `quadratic` and not `quadrature` elements. Both

``````deg = 2
V = VectorFunctionSpace(mesh, "CG", deg)
``````

and

``````W = VectorElement("CG", mesh.ufl_cell(), deg)
V = FunctionSpace(mesh, W)
``````

are equivalent!!

2 Likes

Yeah, because actually I have just got familiar with FEniCS and FEM. I am sorry I did not make my question clear. Thank you for your added up information. This also helps me a lot.