Thank you for comment, @dokken! It was very helpful. I really appreciate it!

*(Apologies for the late reply. )*

In case it helps anyone reading this in the future, here’s my understanding based on the code you answered(bottom of my reply, “MWE”).

To wrap things up, I noticed a difference in the `u`

array results between the code I originally attached and the code you replied to, and between **projection** and **interpolation**.

Despite being a simple mesh(I changed the mesh size 10 \times 10 to 2 \times 2 to make it easier to see).

Code:

```
from dolfin import *
mesh = UnitSquareMesh(2, 2)
CG1 = FunctionSpace(mesh, "CG", 1)
x, y = SpatialCoordinate(mesh)
u = x**2 + y**2
print(f"size of u vector: {project(u, CG1).vector().get_local().size}")
print(f"u vector: {project(u, CG1).vector().get_local()}")
print(f"interpolate u: {interpolate(Expression("x[0]*x[0]+x[1]*x[1]", degree=2), CG1).vector().get_local()}\n")
```

Output:

```
size of u vector: 9
u vector: [ 0.93571429 0.16428571 1.16428571 -0.09285714 0.42142857 1.90714286
0.16428571 1.16428571 0.93571429]
interpolate u: [1. 0.25 1.25 0. 0.5 2. 0.25 1.25 1. ]
```

There was a difference.

(When I say “gradient operation” later, I mean the same thing as in the code I referenced, but only interpolating the `u`

part, or running it with a difference of `u=x**2+y**2`

)

For me, intuitively, the interpolated result of the two values looks right, but when I did a “gradient operation” with `interpolate u`

, it was different from the exact solution.

On the other hand, when I checked the `u`

array through `project`

, it intuitively looked like there was an error, but when I did the “gradient operation”, it matched the exact solution.

In my opinion, defining `u`

as `u=x**2+y**2`

rather than using `interpolate`

seems to be wrong to the eye, but I think it is correct considering the operation, **is there any way to improve this**?

*(I was hoping to be able to use projection and still see results that were reasonably consistent with the interpolation results)*

MWE:

```
from dolfin import *
mesh = UnitSquareMesh(2, 2)
CG1 = FunctionSpace(mesh, "CG", 1)
Vh = VectorFunctionSpace(mesh, "CG", 1)
Th = TensorFunctionSpace(mesh, "CG", 1)
x, y = SpatialCoordinate(mesh)
u = x**2 + y**2
# Hmm... different!
print(f"size of u vector: {project(u, CG1).vector().get_local().size}")
print(f"u vector: {project(u, CG1).vector().get_local()}")
print(f"interpolate u: {interpolate(Expression("x[0]*x[0]+x[1]*x[1]", degree=2), CG1).vector().get_local()}\n")
grad_u = grad(u)
grad_u_ex = as_vector((2*x, 2*y))
print(f"grad_u size: {project(grad_u, Vh).vector().size()}")
print(f"grad_u vector: {project(grad_u, Vh).vector()[:]}")
print(f"x of grad_u size: {project(grad_u.dx(0), Vh).vector().get_local().size}")
print(f"x of grad_u vector: {project(grad_u.dx(0), Vh).vector().get_local()}")
print(f"y of grad_u size: {project(grad_u.dx(1), Vh).vector().get_local().size}")
print(f"y of grad_u vector: {project(grad_u.dx(1), Vh).vector().get_local()}")
grad_error = assemble(inner(grad_u - grad_u_ex, grad_u - grad_u_ex) * dx)
print(f"error between grad_u and grad_u_ex: {grad_error}\n")
Hessian_u = grad(grad_u)
Hessian_u_ex = 2*Identity(2)
print(f"Hessian size: {project(Hessian_u, Th).vector().get_local().size}")
print(f"Hessian vector: {project(Hessian_u, Th).vector().get_local()}")
hes_error = assemble(inner(Hessian_u - Hessian_u_ex, Hessian_u - Hessian_u_ex) * dx)
print(f"error between Hessian_u and Hessian_u_ex: {hes_error}")
```

Output:

```
size of u vector: 9
u vector: [ 0.93571429 0.16428571 1.16428571 -0.09285714 0.42142857 1.90714286
0.16428571 1.16428571 0.93571429]
interpolate u: [1. 0.25 1.25 0. 0.5 2. 0.25 1.25 1. ]
grad_u size: 18
grad_u vector: [-9.58434720e-17 2.00000000e+00 2.51534904e-17 1.00000000e+00
1.00000000e+00 2.00000000e+00 8.32667268e-17 6.44227717e-17
1.00000000e+00 1.00000000e+00 2.00000000e+00 2.00000000e+00
1.00000000e+00 7.13714802e-17 2.00000000e+00 1.00000000e+00
2.00000000e+00 -8.32667268e-17]
x of grad_u size: 18
x of grad_u vector: [2. 0. 2. 0. 2. 0. 2. 0. 2. 0. 2. 0. 2. 0. 2. 0. 2. 0.]
y of grad_u size: 18
y of grad_u vector: [0. 2. 0. 2. 0. 2. 0. 2. 0. 2. 0. 2. 0. 2. 0. 2. 0. 2.]
error between grad_u and grad_u_ex: 0.0
Hessian size: 36
Hessian vector: [2. 0. 0. 2. 2. 0. 0. 2. 2. 0. 0. 2. 2. 0. 0. 2. 2. 0. 0. 2. 2. 0. 0. 2.
2. 0. 0. 2. 2. 0. 0. 2. 2. 0. 0. 2.]
error between Hessian_u and Hessian_u_ex: 0.0
```