solve the Poisson's equation with boundary conditions for the mesh

Question

Hello, it my first time of using Fenics, my issue is: try to apply Poisson's equation with boundary conditions (demo_poisson.py) to my mesh (for instance heat the Cube). I have made Cube (Height=10mm, Length=10mm, Width=10mm) through FreeCad, the mesh was generated by Gmsh. But to solve the problem I have faced obstacles which related to the boundary conditions, at least I think so. I hope you all can help me, and fully explain what I doing wrong. Thank you.

The code is:

    """This demo program solves Poisson's equation

    - div grad u(x, y) = f(x, y)

on the unit square with source f given by

    f(x, y) = 10*exp(-((x - 0.5)^2 + (y - 0.5)^2) / 0.02)

and boundary conditions given by

    u(x, y) = 0        for x = 0 or x = 1
du/dn(x, y) = sin(5*x) for y = 0 or y = 1
"""

# Copyright (C) 2007-2011 Anders Logg
#
# This file is part of DOLFIN.
#
# DOLFIN is free software: you can redistribute it and/or modify
# it under the terms of the GNU Lesser General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# DOLFIN is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Lesser General Public License for more details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
#
# First added:  2007-08-16
# Last changed: 2011-06-28

# Begin demo

from dolfin import *

# Create mesh and define function space 
mesh = Mesh("qwerty.xml")
V = FunctionSpace(mesh, "Lagrange", 1)

# Define Dirichlet boundary (x = 0 or x = 1)
def boundary_1(x):
    return x[2] < 10.0 * DOLFIN_EPS
def boundary_2(x):
    return x[2] > 10.0 - 10.0 * DOLFIN_EPS 

# Define boundary condition
u_1 = Constant(1.0)
bc_1 = DirichletBC(V, u_1, boundary_1)
u_2 = Constant(0.0)
bc_2 = DirichletBC(V, u_2, boundary_2)
bc = [bc_1, bc_2]

# Define variational problem
u = TrialFunction(V)
v = TestFunction(V)
f = Expression("10*exp(-(pow(x[0] - 0.5, 2) + pow(x[1] - 0.5, 2)) / 0.02)")
a = inner(grad(u), grad(v))*dx
L = f*v*dx

# Compute solution
u = Function(V)
solve(a == L, u, bc)

# Save solution in VTK format
file = File("poisson.pvd")
file << u

# Plot solution
plot(u, interactive=True)

Answer 1

Depending on the units of your mesh, you need to set your boundary_1 and boundary_2 functions to return True at the appropriate places. Probably something like:

return abs(x[2] - 10.0) < DOLFIN_EPS 

for the z=10 boundary, for example, and

return x[2] < DOLFIN_EPS

for z=0