Ok, it’s not such a big deal. Thanks a lot man.
Dear @MortezaAram ,
Thanks for sharing your effort! I am also going through the same issue more or less. In my case I am even using just simple linear heat equation to be solved for my simple thermal analysis but still I get same overshoot of temperature near Dirichlet BC when I use some realistic material properties from textbook.
But if I put the value of material properties e.g. rho, cp and thermal_cond as 1 then the overshoot problem disappears all together. Some expert suggested me to use the smoothing function near boundary which was nice and affective solution for my minimal problem but may become complicated for my actual complex problem to implement it overthere. I would like to know how you tackled it as my actual case is near to your issue maybe I am also missing something obvious to be corrected but need another pair of eyes. So far I tried refining the mesh as well as using higher order element but still problem not resolved.
My minimum working problem and snapshot from result is following:
import logging
set_log_level(logging.WARNING)
u_D = Constant(100)
u_initial = Constant(150)
mesh = UnitCubeMesh(15,15,15)
V = FunctionSpace(mesh, 'P', 1)
# Define initial value
u_n = interpolate(u_initial, V)
v = TestFunction(V)
f = Constant(0)
# List
boundary_conditions = {2: {'Dirichlet': u_D}}
class BoundaryZ0(SubDomain):
def inside(self, x, on_boundary):
return on_boundary and near(x[2], 0, 1e-12)
# Mark boundaries
boundary_markers = MeshFunction("size_t", mesh, mesh.topology().dim()-1, 9999)
by0 = BoundaryZ0(); by0.mark(boundary_markers, 2)
# Collect Dirichlet conditions
bcs = []
for i in boundary_conditions:
if 'Dirichlet' in boundary_conditions[i]:
bc_cur = DirichletBC(V, boundary_conditions[i]['Dirichlet'],
boundary_markers, i)
bcs.append(bc_cur)
cp = 500
rho = 4000
thermal_cond = 7
u = Function(V)
utrial = TrialFunction(V)
u.interpolate(u_D)
dt = 1e-2
F = cp*rho*utrial*v*dx + thermal_cond*dt*dot(grad(utrial), grad(v))*dx - cp*rho*(u_n + dt*f)*v*dx
a, L = lhs(F), rhs(F)
# Time-stepping
t = 0
for n in range(200):
t += dt
step = n + 1
solve(a == L, u, bcs)
if n == 0:
# Solver parameters
prm=parameters
info(prm, True)
# Update previous solution
u_n.assign(u)
minima=u.vector().get_local().min()
maxima=u.vector().get_local().max()
print('minima and maxima', minima, maxima, step, t)
xdmff = XDMFFile("thermal_output_mwe/"+format(step, '04')+" .xdmf")
xdmff.write(u)
xdmff.close()
Dear mashsquad777,
Yes, this was quite an issue for me last year and it took me some time to resolve it in the end. But for my case, it came from a simple mistake/ glitch which I’m gonna now share with you.
For the meshing part, I used Iso2Mesh package which renders the nodes unitless. So there was no metric unit for length in my mesh imported to Fenics and this affected drastically the heat transfer equation because some of the parameters like thermal conductivity and heat capacity have the metric unit for length (meter). That’s the reason when we substitute parameters with one the overshoot goes away. The simple solution was for me to scale my mesh according to my real geometry by multiplying the mesh nodes with the appropriate length unit (which was mm in my problem).
I hope this will help you too.
All the best.