# Manufactured solution for poisson type problem **URL:** \ **Category:** dolfinx\ **Created:** [December 12, 2023, 6:21pm UTC](https://fenicsproject.discourse.group/t/manufactured-solution-for-poisson-type-problem/13036 "2023-12-12T18:21:55Z")\ **Posts on this page:** 5\ **Page:** 1
**Author:** ![James\_Dark](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/james_dark/32/2761_2.png) [@James\_Dark](https://fenicsproject.discourse.group/u/James_Dark)\ **Post date:** [December 12, 2023, 6:21pm UTC](https://fenicsproject.discourse.group/t/manufactured-solution-for-poisson-type-problem/13036/1 "2023-12-12T18:21:55Z")
Hello FEniCS experts! I have this script as an MMS test for a Poisson-type problem \nabla(D \cdot \nabla u) + \Gamma = 0 where \Gamma is a source term for an exact solution u = 1 + \sin(2 \pi x). It has been adapted from the Error control tutorial on the [FEniCSx tutorial](https://jsdokken.com/dolfinx-tutorial/chapter4/convergence.html) and a [Previous post](https://fenicsproject.discourse.group/t/manufactured-solution-method-for-multi-material/4059). ```python from mpi4py import MPI from dolfinx import fem import ufl import numpy as np from dolfinx.mesh import create_unit_square from dolfinx.fem.petsc import NonlinearProblem from dolfinx.nls.petsc import NewtonSolver def u_ex(mod): return lambda x: 1 + mod.sin(2 * mod.pi * x[0]) u_ufl = u_ex(ufl) u_numpy = u_ex(np) # mesh my_mesh = create_unit_square(MPI.COMM_WORLD, 50, 50) x = ufl.SpatialCoordinate(my_mesh) elements = ufl.FiniteElement("P", my_mesh.ufl_cell(), 1) V = fem.FunctionSpace(my_mesh, elements) u = fem.Function(V) v = ufl.TestFunction(V) # define bcs u_bc = fem.Function(V) u_bc.interpolate(u_numpy) dofs_L = fem.locate_dofs_geometrical(V, lambda x: np.isclose(x[0], 0)) dofs_R = fem.locate_dofs_geometrical(V, lambda x: np.isclose(x[0], 1)) bc_left = fem.dirichletbc(u_bc, dofs_L) bc_right = fem.dirichletbc(u_bc, dofs_R) bcs = [bc_left, bc_right] D = 1.0 # source term f = ufl.grad(D * ufl.grad(u_ufl(x)))[0, 0] F = 0 F += ufl.dot(D * ufl.grad(u), ufl.grad(v)) * ufl.dx F += f * v * ufl.dx problem = NonlinearProblem(F, u, bcs=bcs) solver = NewtonSolver(MPI.COMM_WORLD, problem) solver.solve(u) ``` I was just wondering if someone would be able to explain the significance of the index `[0, 0]` on the source term `f`. The script only seems to work with it there. Plus I was wondering if there was any way to have this source term more generic i.e is it possible to interpolate it into a function for instance? Thanks in advance, James ---
**Author:** ![dokken](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/dokken/32/1560_2.png) [@dokken](https://fenicsproject.discourse.group/u/dokken)\ **Post date:** [December 12, 2023, 9:40pm UTC](https://fenicsproject.discourse.group/t/manufactured-solution-for-poisson-type-problem/13036/2 "2023-12-12T21:40:09Z")
> [@James\_Dark](#): > > `f = ufl.grad(D * ufl.grad(u_ufl(x)))[0, 0]` This should be `ufl.div(D*….)` not gradient if it is a Poisson problem. > [@James\_Dark](#): > > Plus I was wondering if there was any way to have this source term more generic i.e is it possible to interpolate it into a function for instance? You can interpolate it into a function space, sure, But then you pay the price of approximating the source term in a space it does not lie in (as you are basing your analytical solution on sinus). Is there a specific motivation for wanting to do this? ---
**Author:** ![James\_Dark](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/james_dark/32/2761_2.png) [@James\_Dark](https://fenicsproject.discourse.group/u/James_Dark)\ **Post date:** [December 13, 2023, 10:38am UTC](https://fenicsproject.discourse.group/t/manufactured-solution-for-poisson-type-problem/13036/3 "2023-12-13T10:38:40Z")
Yeah you’re right, should’ve been `ufl.div(D*...)`! Thanks And regarding the interpolation, it’s just for ease of handling the value of the source term within our application, just easier to allow passing of a value of type `fem.Function` rather than allowing a bunch of different ufl objects (`ufl.differentiation.Div` or `ufl.algebra.Product` for instance). But again you’re right it does seem interpolating it onto a function space does seem to introduce approximations, increasing the L2 error. In some testing, I interpolated the value using an `fem.Expression`: ```python f_value = -ufl.div(D * ufl.grad(u_ufl(x))) f_expr = fem.Expression(f_value, V.element.interpolation_points()) f = fem.Function(V) f.interpolate(f_expr) ``` And in doing this it increased the L2 error compared to the analytical solution from 2.5e-08 to 4.7e-08 in this case. Not sure how much this may change with other examples or how acceptable this is… ---
**Author:** ![dokken](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/dokken/32/1560_2.png) [@dokken](https://fenicsproject.discourse.group/u/dokken)\ **Post date:** [December 13, 2023, 10:42am UTC](https://fenicsproject.discourse.group/t/manufactured-solution-for-poisson-type-problem/13036/4 "2023-12-13T10:42:08Z")
As long as you are aware of any errors introduced by the interpolation, I guess you are fine. For analytical problems you can compute the error by computing the L2 error between the interpolated value ad the value evaluated at quadrature points. If you are concerned about type hints, the general type for such “expressions” is `ufl.core.expr.Expr` which can take in any ufl expression. ---
**Author:** ![James\_Dark](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/james_dark/32/2761_2.png) [@James\_Dark](https://fenicsproject.discourse.group/u/James_Dark)\ **Post date:** [December 13, 2023, 2:40pm UTC](https://fenicsproject.discourse.group/t/manufactured-solution-for-poisson-type-problem/13036/5 "2023-12-13T14:40:36Z")
Ok noted, cheers @dokken !