Hi,
I noticed that if i impose a Dirichlet condition of 0 on a domain while solving the diffusion equation, and I set the initial concentration to 0 everywhere, i would get the exact solution of 0 at all times. However, if i repeat the experiment with a nonzero value, I will get small deviations which seem to accumulate over time. Is there any way to tweak the program to avoid this? Below is an example code which reproduces the error.
The reason I'm asking is that I want to simulate a system where there is a small but nonzero concentration source in the system, but I think the noise will disrupt my system.
from dolfin import *
mesh = UnitSquareMesh(11,11)
V = FunctionSpace(mesh, 'CG', 1)
def boundary(x, on_boundary):
return on_boundary
bc_value = 1.0
bc = DirichletBC(V, Constant(bc_value), boundary)
c = interpolate(Constant(bc_value), V)
D = 1.0 # diffusion constant
c_new = TrialFunction(V)
v = TestFunction(V)
dt = 0.0001
form = (c_new - c)*v*dx + dt*inner(D*nabla_grad(c), nabla_grad(v))*dx
(a,L) = system(form)
plot_c = Function(V)
for i in range(1000):
solve(a==L, c)
plot_c.vector()[:] = c.vector().array()-bc_value
plot(plot_c)
