Help in implementing boundary condition of a simple 1D problem.

Question

I have just started learning FEniCS and have used: http://www.scientificpython.net/pyblog/fenics-linear-two-point-bvp to write a script for solving:

           u'' + u = 1

           u(0) = 1, u'(1) = 0

with exact solution,

           u(x) = exp(-x)/ exp(-1) + x + ( -1+exp(-1) )/ exp(-1)

Clearly the weak formulation of the above problem is:

           -(u', v') + (u,v) = (g,v) ; with g = 1

Here is the edited code:

    from dolfin import *

    # definig mesh
    mesh = IntervalMesh(20, 0, 1)

    # definig Function space on this mesh using Lagrange polynoimals of degree 2.
    V = FunctionSpace(mesh, "CG", 2)


    # definign boundary values
     u0 = Constant(0)
     #u0 = Expression("x[0]")

    # this functions checks whether the input x is on the boundary or not.
    def DirichletBoundary(x, on_boundary):
         tol = 1e-14
         return on_boundary and abs(x[0]) < tol

     # Enforcing u = u0 at x = 0
     bc = DirichletBC(V, u0, DirichletBoundary)

     # Setting up the variational problem
     u = TrialFunction(V)
     v = TestFunction(V)
     f = Constant(1)
     g = Constant(1)
     a = -inner(grad(u), grad(v))*dx + inner(u,v)*dx
     L = f*v*dx

     # solving the variational problem.
     u = Function(V)
     solve( a == L, u, bc)

     # plotting solution
     plot(u, interactive = True)

If I use:

   u0 = Constant(1.0),

I get an empty plot.
How can I incorporate the boundary condition u(0) = 1.

Would any body please help?

Edited::::

But If drop the negative sign from the weak form and put f = -1:

                   (u', v') + (u,v) = (f,v) ; with f = -1

then I can run the above code with any Dirichlet boundary condition at x = 0 with the use of:

                  u0 = Constant(2.5), for example.

and this boundary condition is correctly plotted by the FEniCS.

In my knowledge my first stated weak formulation is correct. What am I doing wrong?

Answer 1

Hi, function $u$ that you posted as exact solution is not the exact solution of the problem with right-hand side $f=1.$ In fact $$u`` + u = x - e + 1 + 2e\exp{x}.$$ As far as boundary condition goes, (i) if you require $u(1)=1$ the solution of the problem is $u=1$, the plot might appear empty but that is due to $u$ being constant (check the y-axis)
(ii) with u0 = Expression("x[0]") you are setting the solution to satisfy $u=x|_{x=0}=0$ so you should not expect u(0)=1
(iii) for u(0) = 1 use for example u0 = Expression("1-x[0]")

Comments

Thank you.
I am grateful for your time.