Error setting up the variational form in subdomains 3D

variational formulation
1

Hello, I’m working in composite domain (cylindrical carbon fibre in epoxy matrix) defining the Lamé coefs (lambda and mu) as an expression using c++ for the subdomains, normally I tested my code with a simpler uniform cube for linear elasticity and it works, but after modifying it shows an error when it try to calculate dot(f,v)dx + dot(T,v)ds, the essential part showing the problem:

from fenics import *
from ufl import nabla_grad, nabla_div
import matplotlib.pyplot as plt

f_func = Expression(’-10E4*(pow(x[0], 2) - 25)’, degree=2)

mesh = BoxMesh(Point(-5, -5, 0), Point(5, 5, 1), 80, 80, 10)
V = VectorFunctionSpace(mesh,‘P’, 1)

def boundary_L(x, on_boundary):
tol = 1E-14
return on_boundary and near(x[0], -5, tol)
bc_L = DirichletBC(V, Constant((0, 0, 0)), boundary_L)

def boundary_R(x, on_boundary):
tol = 1E-14
return on_boundary and near(x[0], 5, tol)
bc_R = DirichletBC(V, Constant((0, 0, 0)), boundary_R)

bcs = [bc_L, bc_R]

tol = 1E-14
lamb_0 = 6428571429
lamb_1 = 3.37E11
lambda_ = Expression(‘pow(x[0],2) + pow(x[1],2) <= 2 + tol ? lamb_1 : lamb_0’, degree=0, tol=tol, lamb_1=lamb_1, lamb_0=lamb_0)
mu_0 = 1607142857
mu_1 = 1.444E11
mu_ = Expression(‘pow(x[0],2) + pow(x[1],2) <= 2 + tol ? mu_1 : mu_0’, degree=0, tol=tol, mu_1=mu_1, mu_0=mu_0)

def epsilon(u):
return 0.5*(nabla_grad(u) + nabla_grad(u).T)
def sigma(u):
return lambda_*nabla_div(u)Identity(d) + 2mu_*epsilon(u)

u = TrialFunction(V)
d = u.geometric_dimension()
v = TestFunction(V)
f = (0, f_func, 0)
T = Constant((0, 0, 0))
a = inner(sigma(u), epsilon(v))*dx
L = dot(f, v)*dx + dot(T, v)*ds

u = Function(V)
solve(a == L, u, bcs)

The Error:

Traceback (most recent call last):
File “my_file.py”, line 56, in
L = dot(f, v)*dx + dot(T, v)*ds
File “/usr/lib/python3/dist-packages/ufl/operators.py”, line 158, in dot
a = as_ufl(a)
File “/usr/lib/python3/dist-packages/ufl/constantvalue.py”, line 413, in as_ufl
" to any UFL type." % str(expression))
ufl.log.UFLValueError: Invalid type conversion: (0, Coefficient(FunctionSpace(None, FiniteElement(‘Lagrange’, None, 2)), 0), 0) can not be converted to any UFL type.

I’ve installed FEniCS from source files in my laptop OS Ubuntu 18.04

Thanks by advance for all your comments and help

2

Hi
You can try this:

f_func = Expression('-10E4*(pow(x[1], 2) - 25)', degree=2)
dd =  as_vector([0] + [f_func]+ [0])

T = Constant((0, 0, 0))
a = inner(sigma(u), epsilon(v))*dx
L = dot(dd, v)*dx + dot(T, v)*ds
u = Function(V)

In addition you need to change the linear solver to mumps to avoid RuntimeError:

solve(a == L, u, bcs,solver_parameters={'linear_solver':'mumps'})
1 Like
3

Thanks a lot, that was exactly the problem !!! have a nice day