# How to implement a "project" function correctly?

**URL:** <https://fenicsproject.discourse.group/t/how-to-implement-a-project-function-correctly/5861>\
**Category:** Uncategorized\
**Created:** [May 24, 2021, 3:56pm UTC](https://fenicsproject.discourse.group/t/how-to-implement-a-project-function-correctly/5861 "2021-05-24T15:56:02Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![avishmj](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/avishmj/32/2017_2.png) [@avishmj](https://fenicsproject.discourse.group/u/avishmj)\
**Post date:** [May 24, 2021, 3:56pm UTC](https://fenicsproject.discourse.group/t/how-to-implement-a-project-function-correctly/5861/1 "2021-05-24T15:56:02Z")

</div>

In the following MWE, I am solving the Poisson equation (with the right-hand side multiplied by T) on a disc of radius Sqrt(3). T is a scalar that ranges from 0 to 5. Then in order to apply the solution to every mesh coordinate, I applied the following commands:

```
Y = mesh.coordinates()
X = project(Y,W1)
X += np.vstack(map(u,X))

```

MWE:

```
from fenics import *
from mshr import Circle, generate_mesh
import numpy as np

# Optimization options for the form compiler
parameters["form_compiler"]["cpp_optimize"] = True
ffc_options = {"optimize": True, \
               "eliminate_zeros": True, \
               "precompute_basis_const": True, \
               "precompute_ip_const": True}

class boundary(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary        

class disk(SubDomain):
    def inside(self, x, on_boundary):
        return True

C = Circle(Point(0, 0), sqrt(3))

mesh = generate_mesh(C, 100)
boundary_markers = MeshFunction("size_t", mesh, mesh.topology().dim()-1, 0)
surface_markers = MeshFunction("size_t", mesh, mesh.topology().dim(), 0)

boundary().mark(boundary_markers, 2)
disk().mark(surface_markers, 1)

ds = Measure('ds', domain=mesh, subdomain_data=boundary_markers)
dx = Measure('dx', domain=mesh, subdomain_data=surface_markers)
print(assemble(1*dx(1)), assemble(1*ds(2)))

W1 = FunctionSpace(mesh, "Lagrange", 1)

n = FacetNormal(mesh)
G , mu = 1, 0.1
u_D = Constant(0.0)
d = mesh.geometry().dim()
I = Identity(d)            
T = Constant(0.0)
c = Constant((1.0,0.0))
bc = DirichletBC(W1, u_D, boundary_markers, 2)
print(assemble(1*ds(2)))
for k in range(0,5,1):
# Define variational problem
  u = TrialFunction(W1)
  v = TestFunction(W1)
  f = Constant(-G/mu) # f=-G/mu
  a = dot(grad(u), grad(v))*dx
  L = T*f*v*dx 
  T = T + k
# Compute solution
  u = Function(W1)
  solve(a == L, u, bc)
  Y = mesh.coordinates()
  X = project(Y,W1)
  X += np.vstack(map(u,X))

```

Error:

```
 X = project(Y,W1)
  File "/usr/lib/petsc/lib/python3/dist-packages/dolfin/fem/projection.py", line 129, in project
    L = ufl.inner(w, v) * dx
  File "/usr/lib/python3/dist-packages/ufl/operators.py", line 155, in inner
    b = as_ufl(b)
  File "/usr/lib/python3/dist-packages/ufl/constantvalue.py", line 437, in as_ufl
    raise UFLValueError("Invalid type conversion: %s can not be converted"
ufl.log.UFLValueError: Invalid type conversion: [[-1.73197184 0.01653895]
 [-1.73197184 -0.01653895]
 [-1.73134016 -0.04961083]
 ...
 [0.97718622 1.04771416]
 [0.09944235 0.30137168]
 [0.06897023 0.31762106]] can not be converted to any UFL type.

```

What’s wrong here with my implementation of the “project” function or is there another problem with my code? Any help or suggestion will be helpful.

---

<div class="post-metadata">

**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:** [May 24, 2021, 5:09pm UTC](https://fenicsproject.discourse.group/t/how-to-implement-a-project-function-correctly/5861/2 "2021-05-24T17:09:43Z")

</div>

The project function takes in an `Expression`, `Function`, or `ufl.Form`g, and solves the problem:  
\int\_{\Omega} u \cdot v \mathrm{d}x = \int\_{\Omega} g \cdot v \mathrm{d}x where u is the TrialFunction, v the TestFunction.

What do you mean when you say you want to apply the solution to every mesh coordinate? The field of your solution u is a scalar CG 1 field, while the mesh coordinates is a vector field.

---

<div class="post-metadata">

**Author:** ![avishmj](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/avishmj/32/2017_2.png) [@avishmj](https://fenicsproject.discourse.group/u/avishmj)\
**Post date:** [May 25, 2021, 3:30am UTC](https://fenicsproject.discourse.group/t/how-to-implement-a-project-function-correctly/5861/3 "2021-05-25T03:30:53Z")

</div>

Sorry for the late response.

Yes, you are right. I made a mistake. There should be a Vector function space instead of a function space for the solution u. I have made the changes to my MWE

By “apply the solution to every mesh coordinate”, I meant to say that if u=(u\_1, u\_2) is the solution and suppose (x\_1,y\_1) is one my mesh coordinate then the modified mesh coordinate will be X +=(x\_1 + u\_1,x\_2 + u\_2).

Here is the corrected MWE:

```
from fenics import *
from mshr import Circle, generate_mesh
import numpy as np

# Optimization options for the form compiler
parameters["form_compiler"]["cpp_optimize"] = True
ffc_options = {"optimize": True, \
               "eliminate_zeros": True, \
               "precompute_basis_const": True, \
               "precompute_ip_const": True}

class boundary(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary        

class disk(SubDomain):
    def inside(self, x, on_boundary):
        return True

C = Circle(Point(0, 0), sqrt(3))

mesh = generate_mesh(C, 100)
boundary_markers = MeshFunction("size_t", mesh, mesh.topology().dim()-1, 0)
surface_markers = MeshFunction("size_t", mesh, mesh.topology().dim(), 0)

boundary().mark(boundary_markers, 2)
disk().mark(surface_markers, 1)

ds = Measure('ds', domain=mesh, subdomain_data=boundary_markers)
dx = Measure('dx', domain=mesh, subdomain_data=surface_markers)
print(assemble(1*dx(1)), assemble(1*ds(2)))

W1 = VectorFunctionSpace(mesh, "Lagrange", 1)

n = FacetNormal(mesh)
G , mu = 1, 0.1
u_D = Constant((0.0,0.0))
d = mesh.geometry().dim()
I = Identity(d)            
T = Constant(0.0)
bc = DirichletBC(W1, u_D, boundary_markers, 2)

print(assemble(1*ds(2)))

for k in range(0,5,1):

# Define variational problem
  u = TrialFunction(W1)
  v = TestFunction(W1)
  f = Constant(-G/mu) # f=-G/mu
  a = inner(grad(u), grad(v))*dx
  L = T*f*inner(Constant((1.0,1.0)),v)*dx 
  T = T + k

# Compute solution
  u = Function(W1)
  solve(a == L, u, bc)
  Y = mesh.coordinates()
  Z = project(Y, W1)
  Y += np.vstack(map(u,Z))

```

Still the error reads:

```
Z = project(Y, W1)
  File "/usr/lib/petsc/lib/python3/dist-packages/dolfin/fem/projection.py", line 129, in project
    L = ufl.inner(w, v) * dx
  File "/usr/lib/python3/dist-packages/ufl/operators.py", line 155, in inner
    b = as_ufl(b)
  File "/usr/lib/python3/dist-packages/ufl/constantvalue.py", line 437, in as_ufl
    raise UFLValueError("Invalid type conversion: %s can not be converted"
ufl.log.UFLValueError: Invalid type conversion: [[-1.73197184 0.01653895]
 [-1.73197184 -0.01653895]
 [-1.73134016 -0.04961083]
 ...
 [0.97718622 1.04771416]
 [0.09944235 0.30137168]
 [0.06897023 0.31762106]] can not be converted to any UFL type.
```

---

<div class="post-metadata">

**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:** [May 25, 2021, 5:04am UTC](https://fenicsproject.discourse.group/t/how-to-implement-a-project-function-correctly/5861/4 "2021-05-25T05:04:51Z")

</div>

Use `ALE.move(mesh, u)` as explained in [How to translate a 2D mesh by a vector? - #6 by dokken](https://fenicsproject.discourse.group/t/how-to-translate-a-2d-mesh-by-a-vector/5704/6)