### boundary and initial condition

308
views
0
10 months ago by
I am solving 2D heat equation with Dirichlet Boundary on on side and homogeneous Neumann condition on the other sides. Everything is OK. However, I cannot understand the vtk visualization. In all boundary it shows Dirichlet condition!
from fenics import *
import time
T = 2.0
num_steps = 20
dt = T / num_steps
nx = ny = 30
mesh = RectangleMesh(Point(-2, -2), Point(2, 2), nx, ny)
V = FunctionSpace(mesh, 'P', 1)
tol=1e-5
bot = 'near(x[0],-2, tol)'
def boundary_D(x, on_boundary):
if on_boundary:
if bot:
return True
else:
return False
else:
return False

#Dirichlet condition
u_D = Constant(100.0)
bc = DirichletBC(V, Constant(u_D), boundary_D)
#Neuman condition
g = Constant(0.0)
g= interpolate(g,V)

u_0 = Constant(10)
u_n = interpolate(u_0, V)

# Define variational problem
u = TrialFunction(V)
v = TestFunction(V)
f = Constant(10)
g = Constant(0.9)
F = u*v*dx + dt*dot(grad(u), grad(v))*dx - (u_n + dt*f)*v*dx - (dt)*g*v*ds
a, L = lhs(F), rhs(F)
vtkfile = File('heat_gaussian/solution.pvd')

u = Function(V)
t = 0
for n in range(num_steps):
t += dt
solve(a == L, u, bc)
vtkfile << (u, t)
u_n.assign(u)​

paraview visualization for time 0:
except for bot the other boundaries should not be equal to u0=10??

Community: FEniCS Project

4
10 months ago by
bot is a string, it will be always evaluated as True. Either use it directly in the definition
def boundary_D(x, on_boundary):
if on_boundary:
if near(x[0], -2.0, tol):
return True
else:
return False
else:
return False​

or use more compact implementation

def boundary_D(x, on_boundary):
return on_boundary and near(x[0], -2.0, tol)

or omit the definition at all

bc = DirichletBC(V, u_D, 'near(x[0], -2.0)')
def boundary_D(x, on_boundary):
return on_boundary and near(x[0], -2.0, tol)​