Multi-physics: Thermo-elasticity unphysical answer


213
views
1
3 months ago by
I am trying to solve a static thermo-elasticity. In my case, the displacement field is dependent on temperature but the temperature field is not dependent on displacements. I am fixing the left edge of an unit square domain and applying unit displacement on the right edge. The solutions to codes are converging but not giving the correct meaningful physical results. Please have a look at the code: 
from dolfin import *
from fenics import *

E=1;nu=0.3
mu = E/(2*(1+nu)) # Lame's parameters
lamb = E*nu/((1+nu)*(1-2*nu))
alpha=1 # Thermal coefficient of expansion
beta=alpha*(2*mu+3*lamb)
initial_temp=Constant(0.0)

def clamp_boundary(x,on_boundary):
	return on_boundary and near(x[0],0.0)

def right_boundary(x,on_boundary):
	return on_boundary and near(x[0],1.0)

mesh=UnitSquareMesh(50,50)
plot(mesh)
V = VectorFunctionSpace(mesh, "P" , 1) # displacement function space
Q = FunctionSpace(mesh , "P", 1) # temperature function space
W = MixedFunctionSpace([V, Q])

w = Function(W); (u, p) = split(w)
(v, q) = TestFunctions(W)
bc1=DirichletBC(V,Constant((0,0)),clamp_boundary) # Left edge fixed
bc2=DirichletBC(V.sub(0),Constant((1.0)),right_boundary) # Right edge displacement
bc3=DirichletBC(Q,Constant((0)),clamp_boundary) # left edge temperature
bc4=DirichletBC(Q,Constant((100)),right_boundary) # right edge temperature
bc=[bc1,bc2,bc3,bc4]
d=2;

def epsilon(u):
 return 0.5*(nabla_grad(u)+nabla_grad(u).T)

def sigma(u,p):
 return lamb*nabla_div(u)*Identity(d)+2*mu*epsilon(u)-beta*(p-initial_temp)*Identity(d)

a=inner(sigma(u,p),epsilon(v))*dx + inner(grad(p),grad(q))*dx 
F=a
solve(F == 0, w, bc)
(u,p)=split(w)
plot(w[0], title="Disp", interactive=True)
plot(p, title="Temperature", interactive=True)​

 I expect the linear variations in displacement and temperature along the x-axis. But, the displacement and temperature are; 

 

Thank you. 
Fenics Version; 1.6.0

Community: FEniCS Project

4 Answers


5
3 months ago by
hsk  
Problem in applying Boundary conditions use something like this:
bc1=DirichletBC((W.sub(0)),Constant((0,0)),Clamp_boundary) # Left edge fixed
bc2=DirichletBC(W.sub(0).sub(0),Constant((1.0)),Right_boundary) # Right edge displacement
bc3=DirichletBC(W.sub(1),Constant((0)),Clamp_boundary) # left edge temperature
bc4=DirichletBC(W.sub(1),Constant((100)),Right_boundary) # right edge temperature​
1
Exactly, you need to specify the boundary conditions for the space which you use to find the solution, i.e., W:
from fenics import *

E     = 1.0
nu    = 0.3
mu    = Constant(E/(2.0*(1.0+nu))) # Lame's parameters
lamb  = Constant(E*nu/((1.0+nu)*(1.0-2.0*nu)))
alpha = 1.0 # Thermal coefficient of expansion
beta  = Constant(alpha*(2.0*mu + 3.0*lamb))
init_temp = Constant(0.0)

left_bndr = 'on_boundary && near(x[0], 0.0)'
rght_bndr = 'on_boundary && near(x[0], 1.0)'

mesh = UnitSquareMesh(50, 50)

V = VectorFunctionSpace(mesh, 'P', 1) # displacement function space
Q = FunctionSpace(mesh , 'P', 1) # temperature function space
W = MixedFunctionSpace([V, Q])

bc1 = DirichletBC(W.sub(0), Constant((0.0, 0.0)), left_bndr) # Left edge fixed
bc2 = DirichletBC(W.sub(0).sub(0), Constant(1.0), rght_bndr) # Right edge displacement
bc3 = DirichletBC(W.sub(1), Constant(0.0), left_bndr) # left edge temperature
bc4 = DirichletBC(W.sub(1), Constant(100.0), rght_bndr) # right edge temperature
bcs = [bc1,bc2,bc3,bc4]

w = Function(W)
(u, p) = split(w)
(v, q) = TestFunctions(W)

I = Identity(2)

def epsilon(u): return sym(grad(u))

def sigma(u, p):
 return lamb*div(u)*I + 2*mu*epsilon(u) - beta*(p-init_temp)*I

F = inner(sigma(u, p), epsilon(v))*dx + inner(grad(p), grad(q))*dx

solve(F == 0, w, bcs)
(u, p) = w.split()

plot(u[0], title='Disp -- x-comp')
plot(u[1], title='Disp -- y-comp')
plot(p, title='Temperature')
interactive()
written 3 months ago by Adam Janecka  
1
12 weeks ago by
Thank you hsk and Adam Janecka. It works now. 
0
12 weeks ago by
Gideon  
@hsk, @Adam, @Peter, please what are the best (and necessary) books that a new person who wants to learn numerical analysis, FEA, CFD, need to have and digest in order to get to become proficient and effective.
What are the learning paths and step by step procedures to follow in order to become an expert in writing codes to solve problems?
Hi, I have started FEA with this book, probably this will help. book
written 12 weeks ago by hsk  
Thanks @hsk, which did you use to understand how to write code and other materials on FEA?
written 12 weeks ago by Gideon  
0
12 weeks ago by
Gideon  
@hsk, @Adam, @Peter, please what are the best (and necessary) books that a new person who wants to learn numerical analysis, FEA, CFD, need to have and digest in order to get to become proficient and effective.
What are the learning path and step by step procedures to follow in order to become an expert in writing codes to solve problems?
Please login to add an answer/comment or follow this question.

Similar posts:
Search »